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Abstract 

The  authors  investigate  the  hypothesis  that  phytoplankton  blooms  can  be 
modeled  by  treating  phytoplankton  as  discrete  particles  capable  of  self- 
induced  transport  via  buoyancy  regulation  or  other  behaviors.  A  particle¬ 
tracking  model  is  inserted  into  the  CE-QUAL-ICM  eutrophication  model. 
Phytoplankton  are  quantified  as  carbonaceous  biomass  attached  to  the 
particles.  Kinetics  are  adapted  from  CE-QUAL-ICM.  The  new  model  is 
operated  along  with  existing  hydrodynamic  and  watershed  models  of  the 
Potomac  River  estuary.  Initial  application  is  to  the  winter-spring  1994 
period.  The  new  model  is  compared  to  observations  and  to  a  conventional 
model  of  the  spring  diatom  bloom.  The  particle-based  model  successfully 
computes  a  subsurface  chlorophyll  maximum.  The  model  characteristically 
displays  heterogeneous  spatial  distribution  of  chlorophyll  with  high- 
amplitude  oscillations  at  the  semi-diurnal  period.  The  characteristics  and 
applicability  of  the  particle  approach  are  now  established.  The  model 
requires  optimization  of  parameter  values  to  improve  representation  of  the 
observed  bloom  in  the  Potomac. 
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1  Introduction 

Reliable  prediction  of  phytoplankton  transport  and  production  is  central  to 
the  understanding  and  remediation  of  a  host  of  environmental  problems, 
including  eutrophication  and  the  occurrence  of  harmful  algal  blooms. 
Quantitative  phytoplankton  models  have  existed  for  decades  (Riley  1946), 
yet  virtually  all  examples  are  based  on  the  original  framework  wherein 
phytoplankton  are  a  dissolved  substance  transported  passively  and 
exclusively  by  hydrodynamic  processes.  The  authors  will  test  the  hypothesis 
that  phytoplankton  dynamics  (particularly  the  occurrence  of  blooms)  can  be 
more  accurately  predicted  by  treating  phytoplankton  as  discrete  particles 
capable  of  self-induced  transport  via  buoyancy  regulation  or  other 
behaviors. 

One  of  the  earliest  quantitative  phytoplankton  models  was  presented  by 
Riley  (1946).  Riley’s  model  described  the  seasonal  trends  in  phytoplankton 
biomass  as  the  sum  of  three  processes:  photosynthesis,  respiration,  and 
grazing.  These  processes  were  expressed  as  functions  of  fundamental 
variables,  including:  irradiance,  light  attenuation,  temperature,  mixed- 
layer  depth,  herbivore  abundance,  and  nutrient  availability.  Later  models 
(DiToro  et  al.  1971;  Steele  1974)  significantly  expanded  on  Riley’s 
formulations  while  retaining  the  basic  processes.  These  landmark  studies 
all  described  time  variability  in  phytoplankton  biomass  that  was  at  a  single 
location. 

Spatial  variability  was  added  to  phytoplankton  models  through  division  of 
large  systems  into  well-mixed  boxes  (Kremer  and  Nixon  1978;  Thomann 
and  Fitzpatrick  1982).  The  basic  process  relationships  were  applied  within 
each  box  and  exchange  of  material  between  the  boxes  was  quantified  by 
means  of  long-term  flows  and/or  exchange  coefficients.  The  final  step  in  the 
development  of  modern  phytoplankton  models  was  the  coupling  of  phyto¬ 
plankton  dynamics  to  physics-based  multi-dimensional  hydrodynamic 
models  (e.g.,  Cerco  and  Cole  1993;  Fennel  and  Neuman  1996;  Moll  and 
Radach  2003).  These  state-of-the-art  models  fundamentally  assume  that 
1)  phytoplankton  are  well-mixed  throughout  model  cells,  as  defined  by  the 
computational  grid;  and  2)  phytoplankton  move  passively  with  prevailing 
currents  and  turbulence.  At  present,  computer  codes  for  coupled 
hydrodynamic/phytoplankton  simulations  are  widely  distributed 
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(http://www.mvroms.org:  http://www.epa.gov/athens/wwqtsc/html/efdc.html)  and  large-scale, 
multi-dimensional  simulations  are  commonplace.  Despite  the  spread  of 
technology  and  increase  in  computational  complexity,  however,  the  basics 
of  phytoplankton  models  remain  largely  as  expressed  by  Riley  more  than 
60  years  ago. 

The  current,  advanced  models  still  fail  to  describe  all  observed  phyto¬ 
plankton  processes.  For  example,  the  Potomac  River  Estuary,  a  tributary  of 
the  Chesapeake  Bay  system,  regularly  exhibits  a  spring  diatom  bloom  not 
captured  by  conventional  models  due  to  a  fundamental  departure  from 
current  assumptions  of  phytoplankton  dynamics.  Most  significantly, 
maximum  algal  concentrations  commonly  occur  at  the  bottom  of  the 
estuary,  in  the  total  absence  of  light  (Figure  l).  A  conventional  phyto¬ 
plankton  model  predicts  the  maximum  biomass  should  occur  in  the  photic 
zone,  at  or  near  the  surface.  A  second  quandary  is  that  maximum  surface 
concentration  occurs  mid-estuary,  distant  from  the  oceanic  source  of  the 
diatoms  and  from  the  watershed  nutrient  load  (Figure  2).  The  proposed 
conceptual  model  is  that  the  observed  properties  of  the  bloom  represent  a 
propagation  mechanism  (Figure  3).  Oceanic  diatoms  travel  upstream  with 
prevailing  currents  near  the  bottom.  When  they  reach  the  head  of  the  salt 
intrusion,  they  are  advected  upwards  into  the  light,  forming  a  surface 
bloom.  As  they  are  subsequently  carried  downstream  with  prevailing 
surface  currents,  they  sink  towards  the  bottom  (potentially  as  a  result  of 
internal  buoyancy  control),  forming  a  subsurface  concentration  maximum. 
Near  the  bottom,  they  are  carried  upstream  again,  propagating  the  bloom.  A 
similar  mechanism  was  proposed  by  Tyler  and  Seliger  (1978)  to  explain  the 
occurrence  of  dinoflagellate  blooms  in  upper  Chesapeake  Bay. 

Several  model  innovations  are  necessary  to  test  the  proposed  model.  First, 
phytoplankton  must  be  treated  as  particles  with  the  potential  to  behave,  in 
order  to  maintain  themselves  in  favorable  currents;  this  will  be  achieved  by 
building  upon  existing  particle-tracking  technology.  Second,  the  number  of 
individual  cells  must  be  reduced  to  some  feasible  number,  done  by 
employing  super-individuals  that  represent  a  population  of  individual  cells, 
thus  limiting  the  number  of  particles  to  approximately  100,000.  Last,  the 
particles  must  dynamically  interact  with  the  model  environment;  phyto¬ 
plankton  must  be  able  to  respond  to  potential  cues  from  light,  temperature, 
salinity,  and  other  factors.  The  authors  will  install  the  newly-developed 
algorithms  into  an  existing  comprehensive  eutrophication  model  which  will 
provide  transport  and  ambient  conditions  to  the  modeled  particles. 
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Figure  1.  Elevation  view  of  chlorophyll  along  the  Potomac  River  Estuary,  April  1999.  Distance 
is  indicated  from  the  mouth,  at  left,  to  the  head  of  tide,  at  right.  Note  that  chlorophyll 
concentrations  in  excess  of  60  pg  L 1  occur  at  depths  greater  than  10  m  from  km  60  to  80 
and  that  subsurface  concentration  exceeds  surface  concentration  in  that  region.  Note  also 
the  suggestion  that  chlorophyll  is  borne  by  gravitational  circulation,  from  mouth  to  head  near 
the  bottom  and  from  head  to  mouth  at  the  surface. 


Figure  2.  Aerial  view  of  chlorophyll  concentration  in  the  Chesapeake  Bay  system,  April 
1986.  The  Potomac  River  Estuary  is  in  the  lower,  central  portion  of  the  figure.  Note 
that  the  highest  chlorophyll  concentrations  in  the  entire  system  occur  at  an  upwelling 
region,  located  at  the  head  of  the  salt  intrusion  in  the  Potomac  River. 
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Figure  3.  Proposed  mechanism  to  explain  the  Potomac  River  spring 
phytoplankton  bloom  (elevation  view). 
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Exploratory  steps  in  modeling  phytoplankton  as  particles  provide  a  research 
framework  on  which  to  build.  One  approach  examined  the  occurrence  of 
harmful  algal  blooms  by  treating  plankton  as  passive  particles  (Lanerolle  et 
al.  2006).  This  approach  validated  the  concept  of  treating  plankton  as 
particles,  but  omitted  plankton  behavior  and  interactions  with  the  environ¬ 
ment.  An  “agent-based  model”  (ABM),  in  which  phytoplankton  are  treated 
as  motile  particles,  has  been  used  to  describe  vertical  migration  and  cyst 
formation  as  mechanisms  behind  cyanobacteria  blooms  (Hellweger  et  al. 
2008).  A  complete  “Eulerian-Lagrangian  model,”  which  represents  motile 
plankton  particles  in  an  environmental  continuum  has  been  described 
(Woods  2005)  and  at  least  one  application  to  a  prototype  environment  has 
been  conducted  (Grieco  et  al.  2005).  These  studies  indicate  that  modeling 
phytoplankton  as  particles,  passive  or  active,  is  conceptually  possible  and 
that  the  technique  has  promise.  A  tremendous  amount  of  research,  how¬ 
ever,  remains  in  developing  the  requisite  algorithms  to  capture  phyto¬ 
plankton  dynamics,  coupling  with  spatially  and  temporally  resolved  water 
quality  models,  and  investigating  the  applicability  of  these  approaches  to 
various  systems  where  the  current  advanced  models  do  not  sufficiently 
predict  phytoplankton  distribution. 
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2  ICM  Particle  Tracking1 

Introduction 

Particle  tracking  methods  can  generally  be  divided  into  two  categories 
depending  upon  the  way  in  which  they  treat  the  model  grid.  One  approach 
is  to  perform  the  calculations  on  the  physical  plane  of  the  grid,  while  the 
other  is  to  perform  the  calculations  on  a  transformed  numerical  grid.  The 
particle  tracking  application  by  Chapman  et  al.  (1994)  used  a  transformed 
grid  to  perform  all  calculations.  The  calculations  on  the  transformed  grid 
may  appear  simpler  initially;  however,  the  drawbacks  with  this  approach 
are  the  numerical  complexity,  the  need  to  transform  the  grid  both  before 
and  after  calculations,  and  the  complication  of  the  code.  For  the  ICM 
implementation,  the  calculations  are  performed  within  the  physical  plane 
through  the  use  and  application  of  general  algebra,  geometry,  and  physics 
equations. 

Theory  and  Implementation 

Although  the  ICM  grid  is  unstructured,  the  orthogonal  box  approach  can 
still  be  used  with  some  slight  modifications.  Consider  a  particle  at  position  P 
in  the  box  shown  in  Figure  4.  This  figure  is  shown  in  only  the  x-y  plane  for 
simplicity. 

Although  none  of  the  cell  faces  are  orthogonal,  an  orthogonal  box  can  be 
created  using  the  points  at  which  the  particle  would  exit  the  cell  if  moving 
along  the  x  and  y  axes.  This  is  illustrated  in  Figure  5. 

The  dashed  lines,  A-B  and  C-D,  show  the  paths  that  the  particle  would  take 
if  it  moved  strictly  along  either  the  x  or  y  axis.  The  points  at  which  the 
particle  would  exit  the  grid  cell  along  these  paths  are  then  used  to  construct 
an  orthogonal  box  as  shown  in  Figure  6. 

Since  the  box  shown  in  Figure  6  is  perfectly  orthogonal  and  the  movement 
of  the  particle  along  each  axis  is  performed  separately,  the  velocities  at 
points  A,  B,  C,  and  D  can  be  used  along  with  the  formulas  from  Zheng  and 
Bennett  (2006a,  2006b)  to  determine  the  particle’s  final  location  through 


1  Portions  of  this  chapter  are  adapted  from  material  provided  by  Mr.  Scott  Fant,  formerly  of  SpecPro  Inc. 
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Figure  4.  Initial  position  of  a  single  particle 
in  an  unstructured  grid  cell. 


Figure  5.  Particle  paths  along  x  and  y  axis. 
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Figure  6.  Orthogonal  box  for  bounding 
particle  movement. 


integration.  Although  the  lectures  from  Zheng  and  Bennett  (2006a,  2006b) 
show  the  use  of  these  methods  for  tracking  particles  in  groundwater  flow, 
the  same  equations  should  apply  to  movement  in  surface  water. 

The  velocities  at  points  A,  B,  C,  and  D  are  determined  by  finding  the  vector 
components  of  the  velocities  through  the  faces.  Since  the  flows  are  given 
perpendicular  to  the  flow  faces,  the  angle  of  each  face  relative  to  its 
appropriate  axis  is  required  to  determine  the  components  of  the  vector 
along  the  x  and  y  axes.  The  points  at  the  cell  corners  are  used  to  determine 
the  equations  of  the  facial  lines  in  the  x-y  plane.  This  gives  the  slope, 
y-intercept,  and  facial  angle  for  each  flow  face.  The  equations  for  the  facial 
lines  are  also  used  to  determine  whether  a  particle  has  passed  through  one 
of  the  flow  faces. 


Once  the  velocities  are  known,  the  following  equation  from  Zheng  and 
Bennett  (2006a,  2006b)  is  used  to  determine  the  amount  of  time  required 
for  the  particle  to  cross  the  cell  boundary  in  the  x-direction.  Similar 
equations  are  used  for  the  y  and  z  axes: 


Afx  Aln 


V 


x2 


V 


xl 


(1) 
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in  which: 

vXp  =  Particle  velocity  in  the  x-direction  at  current  location 
Vxi  =  Particle  velocity  in  the  x-direction  at  point  A  of  the  bounding 
box 

Vx2  =  Particle  velocity  in  the  x-direction  at  point  B  of  the  bounding 
box 


Ax  =  Distance  between  points  A  and  B  of  the  bounding  box 

This  assumes  that  the  flow  will  be  in  the  positive  x-direction.  If  the  flow  is 
in  the  negative  x-direction,  all  occurrences  of  uXi  in  Equation  l  should  be 
replaced  with  vX2,  and  vice-versa. 

The  time  that  is  calculated  in  Equation  l  is  used  to  determine  the  maximum 
time  step.  If  the  times  calculated  for  flow  along  all  axes  are  more  than  the 
current  model  time  step,  the  model  time  step  is  used  in  the  calculations.  If 
one  or  more  of  the  calculated  values  is  less  than  the  model  time  step,  the 
lowest  value  is  used  in  the  calculations.  After  movement  of  the  particle,  a 
loop  is  performed  and  calculations  are  made  for  the  remainder  of  the  model 
time  step. 

The  equation  used  to  determine  the  particle’s  final  location  is  given  below. 

x(t2)  =  Xi  +  —  [vxp  (t± )exp  (AxAt)  -  vxl ]  (2) 


in  which: 

x(t2)  =  Particle  location  at  the  end  of  the  time  step 
Xi  =  Particle  location  at  the  beginning  of  the  time  step 
At  =  Particle  tracking  time  step 

After  performing  this  calculation  for  movement  along  each  axis,  the 
particle’s  final  location  at  the  end  of  the  time  step  will  be  determined. 
However,  due  to  the  irregular  shape  of  the  grid  cells,  it  is  possible  for  a 
particle  to  cross  a  grid  cell  boundary  without  crossing  any  of  the  bounding 
box  faces.  The  equations  for  the  facial  lines  are  used  to  determine  whether  a 
particle  is  outside  of  the  grid  cell.  If  the  particle  crosses  a  solid  boundary, 
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the  particle  position  is  reset  to  the  edge  of  the  boundary.  If  the  particle 
crosses  into  a  new  grid  cell,  the  containing  cell  information  is  updated.  If  the 
particle  moves  out  of  the  grid  through  a  bounding  flow  face,  the  containing 
cell  is  set  to  “o”  and  calculations  for  the  particle  will  no  longer  be 
performed. 

Diffusion  and  Dispersion 

The  effects  of  turbulent  diffusion  and  dispersion  are  translated  into 
random  displacement  in  the  particle  model.  The  random  displacement  is 
added  to  the  deterministic  location  derived  from  advection.  For  spatially- 
uniform  diffusion,  the  displacement  is  (Hellweger  and  Bucci  2009): 

A x  =  R.42KAt  (3) 


in  which: 

Ax  =  Particle  displacement 

R  =  Random  number  from  a  standard  normal  distribution 
K  =  Diffusion  or  dispersion 
At  =  Particle  tracking  time  step 

Axial  and  transverse  dispersion  are  modeled  with  a  constant  value,  12  m2  s_1, 
and  application  of  Equation  3  is  straightforward.  Vertical  turbulent  diffu¬ 
sion  is  spatially  varying  and  influenced  by  turbulent  shear,  mixing  length, 
and  density  stratification.  More  sophisticated  formulae  for  addressing 
spatially  varying  diffusion  are  available  (Ross  and  Sharpies  2004)  although 
their  implementation  can  be  complex  on  the  unstructured  ICM  computa¬ 
tional  grid.  As  a  first  approximation,  the  formula  for  spatially  uniform 
diffusion  is  employed  using  a  characteristic  value  K  =  lcrs  m2  s_1.  Investiga¬ 
tion  of  the  effects  of  spatially  varying  diffusion  is  left  for  additional 
investigation  following  completion  of  the  initial  stage  of  this  investigation. 

Caveats 

Although  the  basic  implementation  may  appear  simple  in  nature,  there  are 
several  issues  that  must  be  considered  and  dealt  with  as  appropriate.  When 
a  particle  is  near  the  edge  of  a  cell,  the  bounding  box  may  be  rather  small. 
This  could  cause  the  particle  to  approach  the  boundary  asymptotically,  in 
which  case  the  particle-tracking  module  may  appear  to  loop  endlessly.  In 
order  to  prevent  this,  the  bounding  box  is  generally  extended  by  a  distance 
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of  two  meters  in  each  direction  to  allow  for  movement  of  the  particle  across 
the  cell  boundary.  Flow  toward  a  solid  boundary  can  also  cause  the  same 
type  of  problem  since  a  particle  that  passes  a  solid  boundary  is  reset  to  the 
edge  of  the  boundary.  To  avoid  this  situation,  the  velocity  in  the  direction  of 
a  solid  boundary  is  set  to  a  value  of  “o”  if  the  particle  is  within  one  meter  of 
the  boundary. 

Another  issue  that  can  cause  problems  is  the  orientation  of  the  grid  cells. 
Due  to  the  unstructured  nature  of  the  grid,  positive  flow  in  the  x-direction 
is  not  always  along  the  x-axis.  This  is  especially  true  in  tributaries  that  may 
curve  dramatically.  This  is  illustrated  below  in  Figure  7. 


Figure  7.  Grid  cells  in  a  tributary  indicating 
flow  in  the  positive  x-direction. 


As  seen  in  this  figure,  all  of  the  grid  cells  are  linked  end-to-end.  However, 
the  direction  of  flow  through  these  faces  changes  due  to  the  irregularities  of 
the  cell  shapes.  This  does  not  impact  transport  of  contaminants  between 
cells,  but  requires  a  routine  to  reorient  the  faces  and  facial  flows  of  the  cell 
in  order  to  perform  the  calculations  accurately  for  spatially  aware  particles. 
This  routine  requires  that  the  grid  cell  coordinates  be  specified  in  a  clock¬ 
wise  manner  in  the  input  file.  Although  it’s  possible  to  write  a  routine  to 
work  with  grid  points  specified  in  the  opposite  direction,  it  is  much  easier  to 
simply  specify  the  data  in  the  required  manner. 
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As  mentioned  previously,  a  check  is  performed  on  the  particle  location  after 
each  time  step  to  determine  if  the  particle  has  passed  through  a  cell  face. 
However,  if  the  particle  passes  through  the  corner  of  the  cell,  the  new 
containing  cell  can  be  difficult  to  determine.  For  this  reason,  another  check 
of  the  particles  location  is  performed  at  the  beginning  of  the  calculations.  If 
the  particle  is  found  to  be  outside  of  the  current  cell,  the  containing  cell  is 
updated,  and  the  test  is  performed  again.  Although  this  works  reasonably 
well,  there  are  certain  cases  for  which  this  may  cause  an  endless  loop.  This 
issue  does  not  occur  frequently  and  has  yet  to  be  resolved.  It  appears  that 
the  calculations  are  not  accurate  enough  to  detect  the  containing  cell  of  a 
particle  that  is  very  close  to  a  facial  boundary.  This  only  seems  to  occur  with 
very  small  grid  cells. 
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3  Phytoplankton  Kinetics 

A  primary  goal  of  this  investigation  is  to  transport  algae  as  particles,  rather 
than  by  analogy  to  a  dissolved  substance.  Aside  from  potential  differences  in 
net  transport,  this  approach  allows  for  incorporation  of  algal  behavior,  such 
as  rising  and  falling  via  buoyancy  regulation.  Algae  are  quantified  as 
biomass  (carbon)  per  particle.  The  number  of  particles  in  the  system  is 
constant.  Changes  in  algal  biomass  are  represented  by  altering  the  biomass 
attached  to  a  particle  rather  than  by  changing  the  number  of  particles.  The 
particle-based  kinetics  are  incorporated  into  the  framework  provided  by  the 
CH3D  hydrodynamic  model  (Johnson  et  al.  1993;  Kim  2012)  and  the  ICM 
eutrophication  model  (Cerco  and  Cole  1993;  Cerco  and  Noel  2012) 

Basic  Equation 

The  basic  equation  for  algal  biomass  is: 

^-  =  (G-Pr).B  (4) 

at 


in  which: 

B  =  Algal  biomass  per  particle  (g  C) 

G  =  Net  growth  rate  (d1) 

Pr  =  Net  predation  rate  (d_1) 

The  computer  code  employs  discrete  arithmetic  so  that  the  solution  to 
Equation  3  is: 


Bt+M  =  £f.(l  +  (G-Pr).Af) 


(5) 


in  which: 

Bt+At  =  Algal  biomass  per  particle  at  time  t+At  (g  C) 

Bl  =  Algal  biomass  per  particle  at  time  t  (g  C) 

At  =  Discrete  model  time  step  (d) 

The  net  growth  rate  is  derived  from  photosynthesis  minus  respiratory 
losses: 
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pB 

G  = - .(1  -  PRSP)  -  BM  (6) 

CChl  v  7 

in  which: 

PB  =  Photosynthetic  rate  (g  C  g_1  Chi  d_1) 

CChl  =  Carbon-to-chlorophyll  ratio  (g  C  g1  chlorophyll  a) 

PRSP  =  Active  respiration  (o  <  PRSP  <  l) 

BM  =  Basal  metabolism 

The  photosynthetic  rate  is  influenced  by  ambient  light,  temperature,  and 
nutrient  concentration  as  described  by  Cerco  (2000),  Cerco  and  Noel 
(2004)  and  Cerco  at  al.  (2010).  These  references  also  describe  the  formula¬ 
tions  for  metabolism,  predation,  and  other  influences  on  algal  biomass. 

From  Particles  to  the  Continuum 

The  rate  of  change  in  biomass  attached  to  a  particle  is  determined  by  the 
environment  surrounding  the  particle;  e.g.,  light  and  temperature.  The  rate 
of  change  also  influences  the  environment;  e.g.,  ambient  nutrient  concen¬ 
tration.  A  means  is  necessary  to  translate  from  the  particle-based  system  to 
the  continuum.  The  approach  employed  here  takes  advantage  of  the  discrete 
computational  grid  employed  by  the  hydrodynamic  model  and  by  the  ICM 
eutrophication  model.  ICM  kinetics  are  based  on  ambient  conditions  and  on 
concentrations  computed  in  individual  computational  cells.  The  cell-based 
conditions  and  concentrations  provide  the  environment  for  the  particles. 
The  particles  are  arrayed  in  a  computational  structure  which  labels  each 
particle  with  its  current  grid  cell.  The  structure  and  labeling  allow  the 
computer  code  to  loop  over  all  particles  or  over  all  cells.  The  code  can  query 
particles:  “What  cell  am  I  in?”  Or  it  can  query  cells:  “What  particles  do  I 
contain?”  The  revised  ICM  code  with  particle-based  kinetics  employs  the 
algal  subroutine  from  the  conventional  ICM  code  to  provide  the  particle 
environment,  compute  growth  and  loss  rates,  and  provide  interactions  with 
quantities  computed  in  the  continuum  such  as  nutrient  concentration.  The 
algal  concentration  is  computed  as: 

Elfi, 
v, 


Cone- 


(7) 
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in  which: 

Conci  =  Algal  concentration  in  grid  cell  i  (g  C  nr3) 
n  =  Number  of  particles  in  grid  cell  i 
Bi,p  =  Algal  biomass  attached  to  particle  P  in  grid  cell  i  (g  C) 

V  =  Volume  of  grid  cell  i 

Particle-based  photosynthetic  rate  is  computed: 

p.B  —  PBm»f(I) 

CID  •/(N) 

CID  •/cn 

CID  (8) 


in  which: 

PiB  =  Photosynthetic  rate  associated  with  particle  i  (g  C  g_1  Chi  d_1) 
PBm  =  Maximum  photosynthetic  rate  under  ideal  conditions  (g  C  g1 
Chi  d1) 

f(I)ciD  =  Effect  of  light  on  photosynthetic  rate  in  cell  CID,  containing 
particle  i  (o  <  f(I)  <  l) 

f(N)ciD  =  Effect  of  nutrients  on  photosynthetic  rate  in  cell  CID, 
containing  particle  i  (o  <  f(N)  <  l) 
f(T)cro  =  Effect  of  temperature  on  photosynthetic  rate  in  cell  CID, 
containing  particle  i  (o  <  f(T)  <  l) 


Environmental  effects  on  particle-based  metabolism  and  similar 
parameters  are  computed  by  analogy  to  Equation  8. 

Kinetics  Rules  for  Particles 

Particle-tracking  models  conventionally  incorporate  rules  for  particle 
behavior.  Particle  position  is  influenced  by  rules  as  well  as  by  hydro¬ 
dynamics.  Rules  are  especially  common  for  determining  particle  behavior 
when  the  trajectory  takes  a  particle  outside  the  model  domain.  The  rules 
associated  with  the  present  particle-tracking  algorithms  are  described  in 
Chapter  2.  Physical  rules  for  particle  behavior  potentially  affect  the  attached 
algal  biomass.  Behavioral  rules  can  also  be  formulated  based  on  influences 
from  the  attached  biomass,  instead  of  physical  principles.  Three  rules  are 
enforced  which  influence  particle  behavior  and  associated  algal  biomass. 
These  rules  are  intended  to  keep  the  number  of  particles  in  the  system 
constant  and  to  ensure  that  each  particle  is  associated  with  viable  algal 
biomass: 
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1.  When  a  particle  leaves  the  system  through  an  open  boundary,  the  particle 
with  the  largest  attached  algal  biomass  is  split  in  two.  Half  the  algal 
biomass  is  assigned  to  each  particle. 

2.  When  a  particle  adheres  to  a  solid  boundary,  the  particle  is  removed  from 
the  system.  The  particle  with  the  largest  attached  algal  biomass  is  split  in 
two.  Half  the  algal  biomass  is  assigned  to  each  particle. 

3.  If  algal  biomass  associated  with  a  particle  declines  to  less  than  1%  of  the 
initial  value,  the  particle  is  removed  from  the  system.  The  particle  with  the 
largest  attached  algal  biomass  is  split  in  two.  Half  the  algal  biomass  is 
assigned  to  each  particle. 
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4  The  Potomac  River  Estuary 

The  Potomac  River  Estuary  (Figure  8)  is  a  major  sub-estuary  of  the  larger 
Chesapeake  Bay  system.  The  estuary  extends  190  km  from  the  junction 
with  Chesapeake  Bay  to  the  head  of  tide  at  Washington  DC.  Mean  tide 
range  near  the  mouth  is  0.38  m  and  is  0.84  m  at  Washington  DC.  The 
estuary  is  a  drowned  river  valley  and  the  saline  portion  is  weakly-stratified 
by  Pritchard’s  classification  (Pritchard  1955).  The  primary  freshwater 
source  is  from  the  29,940  km2  upland  watershed  and  enters  at  the  head  of 
tide.  Lesser  volumes  enter  the  estuary  from  the  adjacent  watershed  below 
the  head  of  tide.  Long-term  mean  runoff  at  Washington  DC  is  339  m3  s1. 
This  long-term  mean  is  subject  to  regular  seasonal  fluctuations  as  well  as 
extremes  due  to  flood  and  drought. 

Salinity  intrudes  roughly  140  km  from  the  mouth  to  the  vicinity  of  Station 
TF2.4  (Figure  8).  The  EPA  Chesapeake  Bay  Program  (CBP)  divides  the 
estuary  into  three  segments  based  on  characteristic  salinity,  geometry,  and 
other  features.  The  tidal  fresh  segment  corresponds  to  the  region  occupied 
by  the  stations  designated  “TF”  in  Figure  8.  The  oligohaline  segment  is 
occupied  by  the  stations  designated  “RET.”  The  mesohaline  segment  com¬ 
prises  the  remainder  of  the  estuary.  The  spring  algal  bloom  that  is  the 
subject  of  this  investigation  is  confined  to  the  oligohaline  and  mesohaline 
regions.  Physical  characteristics  of  these  regions  are  summarized  in  Table  1. 

The  Chesapeake  Bay  Environmental  Model  Package  (CBEMP) 

The  particle-tracking  model  is  inserted  into  the  existing  CBEMP.  The 
CBEMP  consists  of  three  independent  models:  a  watershed  model  (WSM), 
a  hydrodynamic  model  (HM)  and  a  eutrophication  model  (WQM).  The 
WSM  (Shenk  and  Linker  2012)  provides  distributed  flows  to  the  HM  and 
nutrient  and  solids  loads  to  the  WQM.  The  HM  (Johnson  et  al.  1993;  Kim 
2012)  computes  three-dimensional  intra-tidal  transport  and  supplies 
transport  parameters  to  the  WQM  on  an  hourly  basis.  The  WQM  (Cerco  et 
al.  2010)  computes  algal  biomass,  nutrient  cycling,  and  dissolved  oxygen, 
as  well  as  numerous  additional  constituents  and  processes. 

Both  the  HM  and  WQM  operate  on  a  three-dimensional  grid  that 
encompasses  the  entire  Chesapeake  Bay  system.  For  this  study,  the 
Potomac  River  portion  of  the  grid  (Figure  9)  was  extracted  along  with  the 
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Table  1.  Physical  characteristics  of  Potomac  River  Estuary  segments.  These  are  derived  from 
corresponding  properties  of  the  model  computational  grid. 


Tidal  Fresh  (TF) 

Oligohaline  (RET) 

Mesohaline  (MH) 

Length  (km) 

61.4 

44.0 

77.6 

Surface  Area  (km2) 

144 

207 

837 

Volume  (106  m3) 

700 

1089 

6542 

Maximum  Depth  (m) 

14.4 

12.8 

20.5 
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Figure  9.  The  Potomac  River  portion  of  the  Chesapeake  Bay  computational  grid.  Depths  are 
shown  from  shallow  (red)  to  deep  (blue). 


associated  hydrodynamics.  An  adjacent  portion  of  Chesapeake  Bay  was 
included  so  that  downstream  boundary  conditions  could  be  specified  at  a 
distance  sufficient  to  minimize  influence  on  the  upstream  portions  of  the 
system.  The  resulting  grid  extended  the  190  km  length  of  the  Potomac  River 
estuary  (Figures  10, 11)  and  incorporated  numerous  embayments  and 
tributaries.  The  grid  consisted  of  1739  surface  elements  (650  m  x  1,300  m  x 
1.5  m)  and  8,965  total  elements.  A  seven-year  period,  1994-2000,  was 
simulated  continuously  using  hydrodynamic  time  steps  of  30  seconds  and 
water  quality  time  steps  of  200  seconds. 
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Figure  10.  The  surface  plane  of  the  independent  Potomac  River  computational  grid  showing 

CBP  sample  stations. 
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Potomac  River  Longitudinal  Transect 
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Figure  11.  Elevation  view  of  the  Potomac  River  computational  grid  from  the  mouth  (km  0)  to 

the  head  of  tide  (km  190). 

The  extracted  model  was  checked  against  the  original  model  and  against 
data  to  ensure  that  the  extraction  was  performed  correctly.  This  model 
version  was  used  as  a  starting  point  for  the  particle-tracking  algorithms  and 
associated  algal  kinetics.  The  formulation  and  validation  of  the  WQM  have 
been  described  elsewhere  (Cerco  et  al.  2010;  Cerco  and  Noel  2012)  and  will 
not  be  repeated.  The  authors  concentrate  here  on  new  developments 
associated  with  treatment  of  algae  as  particles. 

Water  Quality  Database 

The  CBP  conducts  a  monitoring  program  throughout  Chesapeake  Bay  and 
tributaries  (Chesapeake  Bay  Program  2012a).  Stations  are  monitored  at 
monthly  intervals  with  measures  conducted  in-situ  (temperature,  salinity, 
dissolved  oxygen,  light  attenuation)  and  samples  collected  at  various 
depths  for  later  analysis  (chlorophyll,  multiple  forms  of  organic  carbon, 
organic  and  inorganic  nitrogen,  organic  and  inorganic  phosphorus,  and 
suspended  solids).  Observations  at  the  stations  designated  “TF,”  “RET,” 
and  “LE”  (Figure  8)  for  the  modeling  period  were  retrieved  from  an  on- 
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line  database  and  subjected  to  various  summary  procedures  for 
comparison  with  the  model. 

Phytoplankton  Database 

The  CBP  also  conducts  a  Living  Resources  Monitoring  Program  which 
includes  a  phytoplankton  component  (Chesapeake  Bay  Program  2012b). 
The  phytoplankton  monitoring  program  is  designed  to  detect  and  monitor 
changes  in  phytoplankton  abundance  and  taxonomic  composition  in 
relation  to  water  quality  conditions  in  Chesapeake  Bay.  Sampling  is 
conducted  at  roughly  monthly  intervals,  although  sampling  intervals  and 
protocols  have  varied  over  the  course  of  the  program.  Relevant  phyto¬ 
plankton  stations  in  the  Potomac  River  are  at  RET2.2  and  LE2.2  (Figure  8). 
Station  CB5.2,  in  the  mainstem  of  the  bay,  provides  information  on  the 
region  encompassing  the  junction  of  the  Potomac  River  with  the  bay. 
Samples  were  collected  above  and  below  pycnocline  at  RET2.2  until  1989, 
after  which  a  single  whole-water  column  sample  was  collected.  Samples 
were  collected  above  and  below  pycnocline  at  LE2.2  and  CB5.2  until  1995 
when  only  a  surface  composite  sample  was  analyzed. 

Phytoplankton  Biomass  and  Species 

Phytoplankton  counts,  by  species,  comprised  the  original  data  contained  in 
the  Phytoplankton  Database.  These  were  converted  to  biomass  as  carbon, 
through  multiplication  by  conversion  factors  (pg  C  cell'1)  associated  with  the 
database.  The  three  diatom  species  with  the  greatest  biomass  during  the 
spring  months  at  the  three  relevant  stations  were  identified  as  Cerataulina 
Pelagica  (CP),  Chaetoceros  Subtilis  (CS),  and  Pseudo-Nitzschia  Pungens 
(PP).  CP  and  PP  were  largely  restricted  to  the  spring  months,  roughly 
February  -  May,  while  CS  was  present  year-round.  The  biomass  of  the  three 
species  in  each  observation  was  summed  to  give  a  total  diatom  biomass  for 
the  spring  bloom  species.  The  total  was  summarized  by  station,  month,  and 
level,  where  appropriate,  for  the  years  1991  -  2000.  The  summary 
(Figure  12)  indicated  the  bloom  occurred  primarily  from  February  -  May. 
While  the  bloom  was  not  prominent,  in  an  average  sense,  in  January,  an 
observation  of  CP  was  sufficient  to  skew  the  median  and  indicated  that 
biomass  of  bloom  magnitude  was  possible  in  that  month.  The  observations 
also  indicated  that  the  bloom  occurred  primarily  downstream  of  RET2.2, 
although  the  upstream  limit  of  the  bloom  could  not  be  determined  due  to 
lack  of  observations  between  RET2.2  and  LE2.2.  Diatoms  were  abundant 
below  the  pycnocline  with  biomass  exceeding  that  above  the  pycnocline  in 
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Figure  12.  Median  monthly  biomass  for  the  sum  of  three  diatoms  at  three  stations,  1991  - 
2000.  (AP  =  above  pycnocline,  BP  =  below  pycnocline,  RET2.2  is  whole-water  sample). 

some  months.  The  degree  to  which  below-pycnocline  biomass  exceeds 
above-pycnocline  maybe  confounded,  however,  by  the  cessation  of  sub- 
pycnocline  sampling  in  1995. 

Biomass  for  individual  species  was  subjected  to  similar  analyses.  At 
RET2.2,  the  highest  bloom  biomass  occurred  during  the  month  of 
February  and  consisted  almost  exclusively  of  CP  and  PP  (Figure  13).  At 
LE2.2,  the  predominant  bloom  species  was  CP,  which  occurred  primarily 
from  February  to  May,  although  the  first  appearance  was  in  December  and 
one  instance  of  enormous  biomass  was  observed  in  January  (Figure  14). 
Median  below-pycnocline  biomass  of  CP  commonly  exceeded  the  above- 
pycnocline  biomass.  The  bloom  of  PP  commenced  a  month  earlier  than 
CP,  in  December,  and  likewise  ended  a  month  earlier,  in  April.  The  timing 
of  the  bloom  at  CB5.2  was  similar  to  LE2.2,  although  at  this  station  PP 
biomass  rivaled  CP  from  December  through  March  (Figure  15).  Median 
below-pycnocline  biomass  of  both  CP  and  PP  frequently  exceeded  the 
above  pycnocline  biomass. 
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Figure  13.  Median  monthly  biomass  for  each  of  three  diatom  species  at  Station  RET2.2, 1991  - 
2000.  (CP  =  Cerataulina  Pelagica,  CS  =  Chaetoceros  Subtilis,  PP  =  Pseudo-Nitzschia  Pungens) 


Figure  14.  Median  monthly  biomass  for  each  of  three  diatom  species  at  Station  LE2.2, 1991 

-  2000. 
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Figure  15.  Median  monthly  biomass  for  each  of  three  diatom  species  at  Station  CB5.2, 1991 

-  2000. 


Chlorophyll  Observations 

Chlorophyll  observations  were  plotted  for  each  survey  from  November 
through  May  for  the  years  1991  -  2000.  Figures  were  produced  as  eleva¬ 
tions,  or  profiles,  along  the  Potomac  axis  from  RET2.2  to  CB5.1  and  were 
superimposed  on  the  model  grid  to  provide  bathymetry  and  spatial  perspec¬ 
tive.  Results  are  presented  for  three  years,  1997  (December  1996  -  May 
199 7),  1998  (December  1997  -  May  1998),  and  1999  (December  1998  -  May 
1999).  These  years  were  selected  because  they  represent,  respectively,  the 
average,  highest,  and  lowest  November  -  May  flows  in  the  ten-year  series. 

Peak  chlorophyll  concentrations  in  December  (Figure  16)  are  in  the  range 
15  -  20  pg  L1  and  are  concentrated  in  the  region  from  RET2.4  to  LE2.3.  By 
January,  concentrations  increase  to  «  30  pg  L1  and  indicate  the  develop¬ 
ment  of  a  subsurface  chlorophyll  maximum  (Figure  17).  By  February 
(Figure  18),  the  maximum  chlorophyll  concentration  approached  100  pg  L1 
in  1997,  the  subsurface  chlorophyll  maximum  was  readily  apparent,  and  the 
bloom  extended  from  RET2.4  to  LE2.3.  By  March  of  the  years  with  low  to 
average  flow,  the  highest  chlorophyll  concentrations  extend  upstream  to 
RET2.2  and  RET2.3,  although  these  concentrations,  5  to  18  pg  L1,  are  lower 
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than  in  other  months  (Figure  19).  In  the  high-flow  year,  chlorophyll 
concentrations  are  of  greater  magnitude  ~  45  pg  Lr1,  occur  between  LE2.2 
and  CB5.1,  and  continue  to  show  a  subsurface  maximum.  By  April,  peak 
chlorophyll  concentrations  have  recovered  to  40  to  50  pg  Lr1,  and  are  in  the 
lower  estuary  from  LE2.2  to  CB5.1  (Figure  20).  By  late  May,  maximum 
chlorophyll  concentrations  are  13  to  30  pg  L1,  and  spatial  patterns  are 
difficult  to  detect  (Figure  21).  The  biomass  plots  suggest  these  chlorophyll 
concentrations  represent  the  last  significant  presence  of  the  CP  and  PP 
bloom  for  the  season  (Figure  14, 15). 


December  of  bloom  season 
from  December  -  May 


Figure  16.  Chlorophyll  concentrations  (ug  L1)  in  December  of  seasons  of  average,  wet,  and 
dry  hydrology.  Samples  are  superimposed  on  an  elevation  view  of  the  computational  grid. 
Station  names  are  shown  above  the  upper  left  panel. 
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Figure  17.  Chlorophyll  concentrations  (ug  L1)  in  January  of  seasons  of  average,  wet,  and  dry 
hydrology.  Samples  are  superimposed  on  an  elevation  view  of  the  computational  grid.  Station 
names  are  shown  above  the  upper  left  panel. 


Date  -  1MN209 


February  of  bloom  season 
from  December  -  May 


Figure  18.  Chlorophyll  concentrations  (ug  L1)  in  February  of  seasons  of  average,  wet,  and  dry 
hydrology.  Samples  are  superimposed  on  an  elevation  view  of  the  computational  grid.  Station 
names  are  shown  above  the  upper  left  panel. 
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Figure  19.  Chlorophyll  concentrations  (ug  L1)  in  late  March  of  seasons  of  average,  wet,  and 
dry  hydrology.  Samples  are  superimposed  on  an  elevation  view  of  the  computational  grid. 
Station  names  are  shown  above  the  upper  left  panel. 


Figure  20.  Chlorophyll  concentrations  (ug  L1)  in  early  April  of  seasons  of  average,  wet,  and  dry 
hydrology.  Samples  are  superimposed  on  an  elevation  view  of  the  computational  grid.  Station 
names  are  shown  above  the  upper  left  panel. 
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Figure  21.  Chlorophyll  concentrations  (ug  L1)  in  late  May  of  seasons  of  average,  wet,  and  dry 
hydrology.  Samples  are  superimposed  on  an  elevation  view  of  the  computational  grid.  Station 
names  are  shown  above  the  upper  left  panel. 
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5  Potomac  River  Salinity  and  Circulation 

Introduction 

The  Potomac  River  computational  grid  and  hydrodynamics  were  extracted 
from  a  complete  Chesapeake  Bay  grid  as  described  in  Chapter  4.  The 
performance  of  the  Potomac  River  portion  of  the  complete  model  was 
examined  and  documented  by  Cerco  et  al.  (2010).  Comparisons  were  made 
of  computed  and  observed  salinity,  of  computed  and  observed  tidal  and 
current  harmonics,  and  of  computed  and  observed  net  circulation.  Here, 
the  authors  present  a  review  of  the  salinity  and  net  circulation  as  com¬ 
puted  on  the  extracted  grid.  Computed  salinity  provides  an  integrated  view 
of  hydrodynamics,  since  it  depends  on  freshwater  runoff,  tidal  flows, 
wind-driven  circulation,  vertical  turbulence,  and  additional  processes.  Net 
circulation  is  a  principal  factor  in  the  hypothesis  regarding  bloom 
maintenance  presented  in  Chapter  1.  The  presentations  are  intended  to 
validate  the  hydrodynamics,  as  computed  on  the  extracted  grid,  and  to 
illustrate  performance  in  regions  and  periods  of  interest  in  the  present 
application. 

Longitudinal  Salinity  Distribution 

The  Chesapeake  Bay  Program  conducts  bay-wide  sampling  surveys  at 
monthly  intervals.  For  comparison  with  the  model  along  the  longitudinal 
axis,  surface  and  bottom  samples  at  each  station  (Figure  8)  are  averaged 
into  seasons  for  each  year.  The  seasons  relevant  to  the  spring  phytoplankton 
bloom  are  winter  (December  -  February)  and  spring  (March  -  May). 

Results  are  presented  here  for  the  spring  season  in  three  years:  1994, 1996, 
and  1999.  These  years  are  emphasized  in  the  complete  bay  simulation  as 
representative  of  average,  wet,  and  dry  conditions,  respectively.  In  the 
Potomac,  the  December  -  May  average  runoff,  as  measured  at  the  fall  line, 
exceeds  the  long  term  average,  503  m3  s_1,  in  both  1994  (886  m3  s_1)  and 
1996  (805  m3  s'1).  The  year  1999  remains  characterized  as  dry  (236  m3  S’1). 

The  longitudinal  comparisons  present  a  broad  summary  of  the  model’s 
ability  to  represent  axial  circulation  and  vertical  stratification.  In  both  of  the 
high-flow  years,  1994  (Figures  22,  23)  and  1996  (Figures  24,  25),  computed 
and  observed  salinity  penetrates  100  km  up  the  estuary  to  the  location  of 
Station  RET2.2  (Figure  8).  Although  the  winter-spring  flows  in  the 
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two  years  differ  by  only  10%,  observed  and  computed  salinity  in  the  lower 
50  km  of  the  estuary  is  depressed  by  «  2  ppt  in  1994  compared  to  1996. 
Observations  at  RET2.2  show  measureable  salinity  in  1999  (Figures  26,  27) 
and  the  location  of  the  salinity  intrusion,  observed  and  computed,  is  10  to 
20  km  upstream  of  the  location  in  the  high-flow  years.  Observed  salinities 
in  the  lower  75  km  of  the  estuary  are  perhaps  5  ppt  higher  than  in  the  high- 
flow  years.  The  model  provides  good  representation  of  the  observed  surface 
salinities  in  this  low-flow  year,  although  the  computed  bottom  salinities 
near  the  mouth  are  low,  indicating  the  computed  stratification  is  less  than 
observed. 

Salinity  Time  Series 

Salinity  time  series  are  presented  for  three  stations  (Figure  8),  one  at  the 
head  of  salt  intrusion  (RET2.2,  Figures  28,  29),  a  second  station  in  the 
transition  region  (RET2.4,  Figures  30, 31)  and  one  in  the  lower  estuary 
(Figures  32  -  34).  Computed  and  observed  salinities  at  the  nominal  head  of 
salt  intrusion  vary  from  zero  to  8  -  10  ppt  over  the  course  of  the  simulation. 
This  station  was  essentially  freshwater  during  the  winter-spring  seasons  of 
1994  and  1996  but  showed  salinity  up  to  *  2  ppt  during  spring  1999.  The 
two  downstream  stations  show  a  progression  of  computed  salinity  from 


POTOMAC  GRID  (Run66) 
Surface  Salinity  Spring  1994 


Figure  22.  Computed  and  observed  surface  salinity  in  the  Potomac  River  for  spring  (March 
-  May)  1994.  Average  observations  are  shown  as  blue  circles  with  range  as  a  vertical  line. 
Computed  average  is  a  red  line  with  range  shown  as  green  dashed  lines. 
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POTOMAC  GRID  (Run66) 
Bottom  Salinity  Spring  1994 


Figure  23.  Computed  and  observed  bottom  salinity  in  the  Potomac  River  for  spring  (March  - 

May)  1994. 


POTOMAC  GRID  (Run66) 
Surface  Salinity  Spring  1996 


Figure  24.  Computed  and  observed  surface  salinity  in  the  Potomac  River  for  spring  (March 

-  May)  1996. 
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POTOMAC  GRID  (Run66) 
Bottom  Salinity  Spring  1996 


Figure  25.  Computed  and  observed  bottom  salinity  in  the  Potomac  River  for  spring  (March  - 

May)  1996. 


POTOMAC  GRID  (Run66) 
Surface  Salinity  Spring  1999 


Figure  26.  Computed  and  observed  surface  salinity  in  the  Potomac  River  for  spring  (March 

-  May)  1999. 
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POTOMAC  GRID  (Run66) 
Bottom  Salinity  Spring  1999 


Figure  27.  Computed  and  observed  bottom  salinity  in  the  Potomac  River  for  spring  (March  - 

May)  1999. 


POTOMAC  GRID  (Run66) 
Salinity  RET2.2  Surface 
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Figure  28.  Time  series  of  observed  (blue  circles)  and  computed  (red  line)  surface  salinity  at 
RET2.2.  The  simulation  runs  from  1994  (year  0)  to  2000  (year  7). 
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POTOMAC  GRID  (Run66) 
Salinity  RET2.2  Bottom 


Years 


Figure  29.  Time  series  of  observed  (blue  circles)  and  computed  (red  line)  bottom  salinity  at 
RET2.2.  The  simulation  runs  from  1994  (year  0)  to  2000  (year  7). 


POTOMAC  GRID  (Run66) 
Salinity  RET2.4  Surface 


Years 


Figure  30.  Time  series  of  observed  (blue  circles)  and  computed  (red  line)  surface  salinity  at 
RET2.4.  The  simulation  runs  from  1994  (year  0)  to  2000  (year  7). 
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POTOMAC  GRID  (Run66) 
Salinity  RET2.4  Bottom 


Figure  31.  Time  series  of  observed  (blue  circles)  and  computed  (red  line)  bottom  salinity  at 
RET2.4.  The  simulation  runs  from  1994  (year  0)  to  2000  (year  7). 
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Figure  32.  Time  series  of  observed  (blue  circles)  and  computed  (red  line)  surface  salinity  at 
LE2.2.  The  simulation  runs  from  1994  (year  0)  to  2000  (year  7). 
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POTOMAC  GRID  (Run66) 
Salinity  LE2.2  Mid-Depth 


Figure  33.  Time  series  of  observed  (blue  circles)  and  computed  (red  line)  mid-depth  salinity  at 
LE2.2.  The  simulation  runs  from  1994  (year  0)  to  2000  (year  7). 


POTOMAC  GRID  (Run66) 
Salinity  LE2.2  Bottom 


Figure  34.  Time  series  of  observed  (blue  circles)  and  computed  (red  line)  bottom  salinity  at 
RET2.4.  The  simulation  runs  from  1994  (year  0)  to  2000  (year  7). 
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least  in  winter-spring  1994  through  1996,  to  highest  in  winter-spring  1999, 
although  the  model  exaggerates  the  higher  salinity  in  1996  in  comparison  to 
1994  at  Station  RET2.4.  This  station  is  located  *  75  km  from  the  river 
mouth,  in  a  region  of  sharp  longitudinal  salinity  gradient.  A  relatively  small 
spatial  difference  in  computed  and  observed  salinity  translates  into  an 
apparently  large  difference  when  results  are  plotted  at  a  fixed  location. 

Potomac  River  Residual  Currents 

A  one-year  current  meter  mooring,  nearly  coincident  with  the  present 
station  LE2.2  (Figure  8),  was  maintained  and  analyzed  by  Elliott  (1978). 
He  reported  that  the  mean  flow  at  a  3-meter  depth  was  seaward,  while  the 
mean  flows  at  7.6  m  and  12.2  m  were  landward.  Throughout  the  year,  the 
daily-average  flow  at  3  meters  reversed  direction  at  intervals  of  four  or  five 
days.  Daily-average  flows  at  deeper  depths  were  almost  exclusively  land¬ 
ward,  although  occasional  reversals  occurred.  Amplitude  of  the  current 
fluctuations  was  ~  20  cm  s-1.  Model  residual  currents  for  one  year  were 
obtained  at  the  location  and  depths  of  Elliott’s  measures.  Residuals  were 
determined  as  the  25-hour  moving  average  of  the  computed  velocities. 
Results  (Figure  35)  compare  well  with  Elliott’s  findings.  Non-tidal  currents 
near  the  surface  reverse  directions  frequently,  while  currents  at  greater 
depths  are  primarily  landward  but  demonstrate  occasional  flow  reversals. 
Currents  at  all  depths  fluctuate  at  periods  of  five  to  six  days.  The 
magnitude  of  the  fluctuations  is  10  to  20  cm  s1. 
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Figure  35.  Model  residual  currents  in  the  Potomac  River  corresponding  to  Elliott’s  (1978) 

measurements. 
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6  Basic  Particle  Simulation 

Introduction 

The  base  simulation  extends  from  November  1, 1993  through  May  31, 
1994.  The  duration  is  selected  to  encompass  the  spring  bloom  period.  The 
year  1994  corresponds  to  an  average  hydrologic  year  selected  for  analysis 
of  Chesapeake  Bay  model  results,  although  flows  in  the  Potomac  during 
the  interval  are  above  average.  Two  flood  events  in  quick  succession  occur 
30  days  into  the  model  run  (Figure  36).  Another  succession  of  flood 
events,  spaced  two  to  three  weeks  apart,  commences  «  90  days  into  the 
run  and  continues  through  the  termination.  Peak  flows  of  4,000  m3  s1 
occur  circa  Day  152. 


Figure  36.  Fall-line  flow  in  the  Potomac  River  November  1, 1993,  to  May  31, 1994. 


Ten  thousand  particles  are  released  45  days  into  the  run  and  are  spaced 
uniformly  from  top  to  bottom  in  the  channel  of  the  river  between  km  30  and 
80.  Particles  are  distributed  in  model  grid  cells  with  the  number  in  each  cell 
specified  to  yield  4  x  io-6  particles  nr3.  The  release  timing  and  region 
correspond  to  the  typical  appearance  date  and  expanse  of  bloom  algae. 
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Presentation  Formats 

Elevation  Plots 

Elevation  plots  extend  along  a  mid-channel  transect  from  the  mouth  of  the 
river  (km  o)  to  the  head  of  tide  (km  160).  Particles  are  shown  only  when 
they  are  within  the  one-cell  wide  transect  and  are  colored  according  to  the 
mass  of  algae  each  represents.  (The  algal  simulations  are  described  in  a 
subsequent  chapter).  The  elevation  plots  are  originally  rendered  as  anima¬ 
tions  from  which  extractions  are  included  in  this  report.  Figure  37  shows 
computed  particle  distributions  at  four  intervals  ranging  from  the  release  to 
the  end  of  the  simulation.  In  this  simulation,  particles  that  leave  the  model 
domain  through  open  boundaries  are  not  replaced.  Consequently,  the 
number  of  particles  in  the  system  decreases  from  the  initial  release. 

JEW*  45  JDay  6tM2742e*01 


Figure  37.  Elevation  view  of  particle  distribution  at  four  time  intervals:  a)  Day  45  (particle 
release);  b)  Day  60;  c)  Day  152;  d)  Day  212  (end  of  simulation).  The  number  of  particles  in 
the  system  declines  as  particles  exit  across  open  boundaries. 


Plan  Views 

Plan  views  illustrate  the  entire  domain  with  particles  superimposed  upon 
the  computational  grid.  Particles  are  shown  at  all  depths  from  surface  to 
bottom  and  are  colored,  according  to  depth,  to  increase  contrast.  Figure  38 
presents  a  plan  view  corresponding  to  the  elevation  plots  in  Figure  37. 
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Figure  38.  Plan  view  of  particle  distribution  at  four  time  intervals:  a)  Day  45  (particle  release); 
b)  Day  60;  c)  Day  152;  d)  Day  212  (end  of  simulation). 


Tracker  Statistics 

The  tracker  statistics  provide  a  format  for  summarizing  the  particle 
distribution  and  number.  The  number  of  particles  is  summarized,  by  level, 
within  10-km  reaches  of  the  river.  Levels  are  derived  from  conventions 
derived  for  the  mainstem  of  Chesapeake  Bay  and  are: 

•  surface  (Level  I),  o  to  6.7  m  from  surface; 

•  mid-Depth  (Level  II),  6.7  to  12.8  m  from  surface;  and 

•  bottom  (Level  III),  more  than  12.8  m  from  surface. 

The  figures  show  particle  counts,  not  concentrations.  At  the  initiation  of  the 
simulation,  the  preponderance  of  particles  is  in  the  surface  (Figure  39),  even 
though  the  concentration  is  uniform  from  surface  to  bottom.  The  vertical 
particle  distribution  reflects  the  greater  volume  of  water  contained  between 
o  and  6.7  m  depth.  For  this  simulation,  corresponding  to  the  results  shown 
in  Figures  37  and  38,  the  number  of  particles  within  the  model  domain 
declines  from  10,000,  at  release,  to  3,764  at  the  end  of  the  simulation. 
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Figure  39.  Tracker  statistics  at  four  time  intervals:  a)  Day  45  (particle  release);  b)  Day  60; 
c)  Day  150;  d)  Day  210  (end  of  simulation).  The  plots  indicate  the  number  of  particles  at 
10-km  intervals  and  three  depth  levels.  The  Potomac  River  mouth  is  at  km  0.  The  bars  at  this 
location  indicate  the  number  of  particles  between  km  0  and  km  10.  Subsequent  graphics  are 
placed  at  the  lower  end  of  10  km  intervals.  The  number  of  particles  outside  the  river  mouth 
but  inside  the  model  domain  is  plotted  at  -10  km. 


Comparison  to  Dissolved  Substance 

The  transport  of  particles  versus  a  dissolved  substance  was  examined  by 
initiating  a  tracer  release  (salt)  simultaneous  with  the  particle  release.  As 
with  the  particles,  the  tracer  was  distributed  uniformly  from  top-to-bottom 
between  km  30  and  80.  Results  (Figure  40)  indicated  that  less  than  one 
percent  of  the  tracer  remained  in  the  model  domain  90  days  after  the 
release  (Day  135  of  the  simulation).  In  contrast,  65%  of  the  particles 
remained  in  the  system  after  90  days  and  nearly  40%  remained  at  the  end 
of  the  simulation,  167  days  after  the  particle  release. 

Elevation  (Figure  41)  and  plan  views  (Figure  42)  of  salinity  illustrate  the 
rapid  dissipation  of  the  dissolved  substance.  For  comparison  with  these 
results,  particle  distribution  was  converted  to  concentration  through 
summing  the  number  of  particles  in  each  grid  cell  and  dividing  by  cell 
volume  (Figures  43,  44).  Aside  from  the  residence  time,  the  dramatic 
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Figure  40.  The  fraction  of  initial  mass  remaining  following  simultaneous  releases  of 
particles  and  dissolved  substance.  Less  than  1%  of  the  dissolved  substance  remains 
90  days  after  release  versus  65%  of  the  particles. 


Figure  41.  Elevation  view  of  dissolved  substance  concentration  at  four  time  intervals:  a)  Day 
45  (release);  b)  Day  80;  c)  Day  115;  d)  Day  150.  The  substance  has  disappeared  by  Day  150. 
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Figure  42.  Surface  view  of  dissolved  substance  concentration  at  four  time  intervals:  a)  Day  45 
(release);  b)  Day  80;  c)  Day  115;  d)  Day  150. 


Figure  43.  Elevation  view  of  particle  concentration  at  four  time  intervals:  a)  Day  45  (release); 
b)  Day  80;  c)  Day  115;  d)  Day  150.  Note  the  “patchy”  distribution  of  particles  compared  to 

dissolved  substance. 
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Figure  44.  Surface  view  of  particle  concentration  at  four  time  intervals:  a)  Day  45  (release);  b) 

Day  80;  c)  Day  115;  d)  Day  150. 


contrast  between  particles  and  dissolved  substance  is  in  the  “patchy” 
distribution  of  the  particles.  The  patchiness  develops  despite  the  uniform 
initial  distribution  of  the  particles. 

Sensitivity  to  Number  of  Particles  Released 

The  sensitivity  of  results  to  the  number  of  particles  was  examined  by 
doubling  the  initial  release  from  10,000  to  20,000  particles  and  then  by 
doubling  it  again  to  40,000  particles.  The  proportion  of  particles  retained 
and  the  spatial  distribution  of  particles  were  similar  for  all  runs,  as 
evidenced  by  the  tracker  statistics  at  Day  210  (Figure  45).  Increasing  the 
number  of  particles  over  the  initial  10,000  placed  a  slightly  larger  fraction 
in  the  lower  two  layers,  but  the  fraction  of  the  total  was  negligible  in  any 
event.  Particle  retention  was  «  38%  of  the  initial  release  regardless  of  the 
number  released.  Based  on  these  findings,  initial  release  of  10,000  particles 
was  retained  in  all  subsequent  runs. 
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Figure  45.  Sensitivity  to  number  of  particles  released.  Tracker  statistics  are  shown  at  the  end 
of  the  simulation  for  release  of:  a)  10,000  particles;  b)  20,000  particles;  c)  40,000  particles. 
The  fraction  particles  remaining  at  the  end  of  the  simulation  and  their  distribution  are  nearly 

identical  for  all  initial  release  numbers. 


Sensitivity  to  Rules 

A  set  of  kinetics  rules  for  particles  is  described  in  Chapter  3.  These  rules 
are  intended  to  keep  the  number  of  particles  in  the  system  constant  and  to 
ensure  that  each  particle  is  associated  with  viable  algal  biomass.  The  rule 
with  the  greatest  effect  on  particle  distribution  is  the  replacement  of 
particles,  which  leave  through  open  boundaries  (Figure  46).  Under  this 
rule,  the  number  of  particles  in  the  system  remains  constant  at  the  initial 
value.  The  rules  governing  the  replacement  of  particles  that  adhere  to  the 
bottom  or  that  no  longer  represent  viable  algae  diminish  the  number  of 
particles  in  the  system.  With  no  rules  in  force,  *  38%  of  the  particles 
released  are  present  at  the  end  of  the  simulation.  When  particles  adhering 
to  the  bottom  or  containing  no  algae  are  replaced,  particle  retention 
declines  to  *  5%.  This  rule  replaces  an  inert  particle  by  splitting  the 
particle  with  the  greatest  algal  biomass.  The  algorithm  effectively  relocates 
particles  from  sluggish  regions  at  the  bottom  into  the  surface  where  they 
are  subject  to  flushing  from  the  system  by  runoff  events.  Enforcement  of 
all  rules  retains  the  initial  number  of  particles  in  the  system,  as  does 
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JD*  2  I  20554e*02  JDay  2  1205S4e*02 


Figure  46.  Sensitivity  to  Kinetics  Rules  for  Particles  (Elevation  views).  Results  are  presented 
for  particle  distributions  at  the  completion  of  simulations  for  four  rule  sets:  a)  No  rules;  b) 
Particles  which  leave  the  system  are  replaced;  c)  Particles  which  adhere  to  the  bottom  or 
contain  no  viable  algae  are  replaced;  d)  Combination  of  b)  and  c). 

simply  replacing  particles  that  leave  the  system.  With  all  rules  in  place, 
however,  the  axial  particle  distribution  is  more  uniform  and  a  greater 
fraction  of  the  particles  are  below  the  surface,  relative  to  the  run  with 
replacement  only  (Figure  47). 

Base  Run  for  Phytoplankton  Simulation 

The  experiments  with  particle  tracking  indicate  that  particle  transport  is 
fundamentally  different  than  transport  of  a  dissolved  substance.  In 
particular,  the  residence  time  of  particles  released  in  the  center  of  the 
Potomac  Estuary  is  longer  than  the  residence  time  of  a  comparable  release 
of  dissolved  substance  (Figure  40).  Ninety  days  after  the  release,  65%  of 
the  particles  remain  in  the  model  domain  while  the  dissolved  substance 
has  disappeared.  A  second  distinction  is  in  the  “patchy”  distribution  of 
particles  compared  to  the  continuous  distribution  of  dissolved  substance. 
The  particle  experiments,  supplemented  by  initial  experiments  with 
phytoplankton,  help  determine  a  base  particle  simulation  (Figures  48,  49) 
for  subsequent  development  of  the  phytoplankton  model.  The  run  has  the 
following  characteristics: 
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•  Release  10,000  particles  at  Day  45  of  a  212-day  simulation. 

•  Enforce  kinetics  rules  that  govern  the  replacement  of  particles  leaving 
the  system,  adhering  to  the  bottom,  or  representing  no  viable  algae. 

•  Incorporate  particle  dispersion  based  on  modeled  axial  dispersion  and 
characteristic  vertical  turbulent  diffusion. 


Figure  47.  Sensitivity  to  Kinetics  Rules  for  Particles  (Tracker  Statistics).  Results  are  presented 
for  particle  distributions  at  the  completion  of  simulations  for  four  rule  sets:  a)  No  rules; 
b)  Particles  which  leave  the  system  are  replaced;  c)  Particles  which  adhere  to  the  bottom  or 
contain  no  viable  algae  are  replaced;  d)  Combination  of  b)  and  c). 
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JD*  45  JD*  6012742e*0t 


Figure  48.  Elevation  view  of  particle  distribution  for  base  phytoplankton  run.  This  run 
incorporates  the  kinetics  rules  set  and  particle  dispersion.  Results  are  presented  at  four  time 
intervals:  a)  Day  45  (particle  release);  b)  Day  60;  c)  Day  152;  d)  Day  212  (end  of  simulation). 


Figure  49.  Tracker  statistics  for  base  phytoplankton  run.  Results  are  presented  at  four  time 
intervals:  a)  Day  45  (particle  release);  b)  Day  60;  c)  Day  150;  d)  Day  210  (end  of  simulation). 
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7  Algae  as  Particles  -  Spatial  and  Temporal 
Distribution 

Introduction 

The  base  model  run  with  algae  as  particles  accompanies  the  basic  particle 
simulation  described  in  the  previous  chapter.  Ten  thousand  particles  are 
released  in  the  central  portion  of  the  Potomac  River  45  days  after 
commencement  of  a  model  run,  which  extends  from  November  1, 1993, 
through  May  31, 1994.  Each  particle  is  initiated  with  2.5  x  ios  g  algal 
carbon,  yielding  1.0  g  C  m-3  at  initiation  of  the  simulation.  Parameters  for 
the  algal  kinetics  (Table  2)  are  drawn  from  the  spring  diatom  group  of  the 
Chesapeake  Bay  Model  (Cerco  et  al.  2010).  These  values  are  employed  to 
provide  comparison  to  algae  as  simulated  in  conventional  models.  No 
attempt  is  made,  at  this  stage  of  development,  to  optimize  the 
performance  of  the  model  with  algae  as  particles. 


Table  2.  Algal  Kinetics  Parameters 


Symbol 

Definition 

Value 

Units 

ANC 

algal  nitrogen-to-carbon  ratio 

0.167 

g  N  gi  C 

APC 

algal  phosphorus-to-carbon  ratio 

0.0125 

gPgiC 

ASC 

algal  silica-to-carbon  ratio 

0.3 

gSigiC 

BM 

basal  metabolic  rate  of  algae  at  reference  temperature  Tr 

0.01 

d1 

CChl 

algal  carbon-to-chlorophyll  ratio 

75 

gcg-i  Chi 

KHn 

half-saturation  concentration  for  nitrogen  uptake  by  algae 

0.025 

g  N  m-3 

KHp 

half-saturation  concentration  for  phosphorus  uptake  by  algae 

0.0025 

g  P  m-3 

KHs 

half-saturation  concentration  for  silica  uptake  by  algae 

0.03 

g  Si  nr3 

KTb 

effect  of  temperature  on  basal  metabolism  of  algae 

0.032 

oc-i 

KTgl 

effect  of  temperature  below  Tm  on  growth  of  algae 

0.0018 

°C"2 

KTg2 

effect  of  temperature  above  Tm  on  growth  of  algae 

0.006 

°C"2 

Phtl 

predation  rate  on  algae 

0.01 

m3  g-i  C  d1 

PmB 

maximum  photosynthetic  rate 

300 

gCg  i  Chi  d1 

Presp 

photo-respiration  fraction 

0.25 

0  <  Presp  <  1 

Topt 

optimal  temperature  for  growth  of  algae 

16 

°C 

Tr 

reference  temperature  for  metabolism 

20 

°C 

Trpr 

reference  temperature  for  predation 

20 

°C 

Wa 

algal  settling  rate  (conventional  model  only) 

0.5 

m  d1 

a 

initial  slope  of  production  vs.  irradiance  relationship 

8.0  (spring 

gCgi  Chi  (E  m-2)-i 
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Two  conventional  model  runs  were  completed  to  provide  comparisons  with 
the  particle  approach.  The  first  is  a  “Duplicate”  run.  Algae  are  initiated  at  a 
concentration  of  1.0  g  C  m-3  at  the  same  time  and  in  the  same  extent  of  the 
Potomac  River  as  the  particle  release.  Concentration  boundary  conditions  at 
the  open  ends  of  the  model  domain  are  specified  as  zero.  Kinetics  para¬ 
meters  are  identical  to  the  particle  simulation  with  one  exception.  The 
Duplicate  run  employs  an  algal  settling  velocity,  as  in  the  Chesapeake  Bay 
model.  Settling  is  omitted  from  the  base  particle  simulation.  The  Duplicate 
run  is  employed  in  the  majority  of  comparisons  with  the  particle  simulation. 
A  run  referred  to  as  “Original”  was  also  completed.  The  original  run  simu¬ 
lates  spring  diatoms  exactly  as  in  the  Chesapeake  Bay  model.  The  foremost 
characteristic  of  this  run  is  that  diatoms  are  initiated  at  the  open  boundaries 
rather  than  in  the  mid-river.  Advection  and  dispersion  transport  diatoms  to 
regions  that  are  favorable  to  growth. 

Plan  and  Elevation  Views 

Elevation  views  of  algae  as  particles  along  the  axial  transect  are  provided 
in  Figure  48.  The  distribution  of  attached  algae  is  heterogeneous  and 
demonstrates  no  apparent  pattern.  Particles  with  varying  amounts  of 
attached  algae  are  distributed  along  the  transect  and  in  the  vertical 
direction.  Particles  with  widely  different  amounts  of  attached  algae  are 
found  immediately  adjacent  to  each  other.  At  Days  152  and  212,  individual 
particles  with  high  amounts  of  attached  algae  occur  at  depths  greater  than 
20  m.  The  amount  of  algae  attached  to  particles  at  great  depths  can  exceed 
the  amount  of  attached  algae  above,  in  the  photic  zone.  In  this  regard,  the 
simulation  resembles  the  observations  (e.g.,  Figures  18, 19)  in  which 
chlorophyll  concentration  near  the  bottom  exceeds  chlorophyll 
concentration  at  the  surface. 

Views  of  algae  as  particles  in  the  model  surface  layer  (Figure  50)  reflect  the 
properties  apparent  in  the  elevation  views.  The  spatial  distribution  is 
heterogeneous  and  lacks  apparent  pattern.  Circulation  carries  algae  into  the 
smallest  tributaries  where  they  remain  when  the  mainstem  of  the  river  is 
largely  cleared  by  the  flood  event  (Figure  36)  circa  day  152.  At  the  comple¬ 
tion  of  the  simulation,  the  amount  of  algae  attached  to  individual  particles 
varies  through  multiple  orders  of  magnitude  in  the  range  >  104  g  C  particle’1 
but  <  107  g  C  particle1.  The  lower  limit  is  affected  by  the  rule  which  elimi¬ 
nates  particles  with  no  viable  attached  algae;  particles  are  eliminated  when 
attached  algae  falls  below  2,500  g  C  particle1.  Particles  with  more  than  2.5  x 
105  g  attached  algae  demonstrate  growth  above  the  initial  concentration. 
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Figure  50.  Surface  view  of  particle  distribution  and  attached  algal  biomass  for  base 
phytoplankton  run.  Results  are  presented  at  four  time  intervals:  a)  Day  45  (particle  release); 
b)  Day  60;  c)  Day  152;  d)  Day  212  (end  of  simulation). 


The  upper  limit  is  also  influenced  by  the  rules  since  the  particle  with  the 
highest  algal  concentration  is  split  to  compensate  for  particles  which  leave 
the  system  or  are  eliminated. 

Comparison  with  Conventional  Model 

Two  conversions  were  conducted  to  enable  comparison  of  algae  computed 
in  the  particle  model  with  algae  computed  in  the  conventional  model.  The 
two  conversions  were  also  conducted  to  facilitate  comparisons  of  both 
models  with  observations.  First,  attached  algal  carbon  was  converted  to 
volumetric  concentration.  Carbonaceous  biomass  in  each  grid  cell  was 
summed  over  all  the  particles  in  the  cell  and  then  divided  by  cell  volume. 
In  the  second  conversion,  carbon  concentration  in  both  models  was 
converted  to  chlorophyll  concentration  through  division  by  the  model 
carbon-to-chlorophyll  ratio.  Comparisons  were  completed  along  the  mid¬ 
channel  transect  and  in  the  surface  plane. 
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The  elevation  views  of  the  particle  model  (Figure  51)  demonstrate  a 
heterogeneous  spatial  distribution  which  manifests  immediately  upon 
release  of  the  particles.  Chlorophyll  concentrations  near  the  bottom 
commonly  exceed  concentrations  at  lesser  depths  above.  The  surface  layer, 
which  receives  the  greatest  illumination,  is  at  times  devoid  of  algae. 
Chlorophyll  concentrations  computed  in  the  conventional,  Duplicate, 
model  are  much  more  smoothly  distributed  (Figure  52).  On  Day  60,  soon 
after  bloom  initiation,  there  is  an  indication  of  a  subsurface  chlorophyll 
maximum  caused  by  upstream  advection  of  the  initial  uniform  chlorophyll 
distribution.  Otherwise,  the  vertical  chlorophyll  distribution  declines 
monotonically  from  the  illuminated  surface  waters  to  the  bottom.  No 
surface  cells  within  the  bloom  region  are  devoid  of  algae. 

Surface  chlorophyll  concentrations  computed  by  the  particle  model  follow 
the  heterogeneous  pattern  previously  observed  in  the  computed  particle 
distribution  (Figure  53).  As  noted  in  the  particle  distributions,  chlorophyll 
penetrates  into  the  smallest  tributaries  and,  by  the  end  of  the  simulation, 
the  highest  concentrations  are  in  the  tributaries  rather  than  in  the  surface 
waters  of  the  mainstem  Potomac.  The  conventional  model  provides  more 
uniform  surface  chlorophyll  distributions  (Figure  54).  At  the  end  of  the 
simulation,  chlorophyll  concentrations  in  the  mainstem  are  as  high,  or 
higher,  than  in  the  tributaries.  The  contrast  between  the  uniform  surface 
distribution  in  the  conventional  model  and  the  patchy  distribution  in  the 
particle  model  is  striking. 
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Figure  51.  Elevation  view  of  chlorophyll  concentration  derived  from  the  model  of  algae  as 
particles.  Results  are  shown  at  four  time  intervals:  a)  Day  45  (release);  b)  Day  60;  c)  Day  152; 
d)  Day  212  (end  of  simulation).  Note  the  “patchy”  distribution  of  chlorophyll  and  the 
occurrence  of  high  concentrations  near  the  bottom  of  the  river. 


Figure  52.  Elevation  view  of  chlorophyll  concentration  computed  by  the  conventional  algal  model. 
Results  are  shown  at  four  time  intervals:  a)  Day  45  (release);  b)  Day  60;  c)  Day  152;  d)  Day  212 
(end  of  simulation).  Note  the  monotonic  decline  of  chlorophyll  from  surface  to  bottom. 


ERDC/EL  TR-13-13 


55 


Figure  53.  Surface  view  of  chlorophyll  concentration  derived  from  the  model  of  algae  as  particles 
(Note  the  use  of  a  log  scale).  Results  are  shown  at  four  time  intervals:  a)  Day  45  (release);  b)  Day 
60;  c)  Day  152;  d)  Day  212  (end  of  simulation).  At  the  end  of  the  simulation,  the  highest 
chlorophyll  concentrations  are  in  small  tributaries. 


Figure  54.  Surface  view  of  chlorophyll  concentration  computed  by  the  conventional  model. 
Results  are  shown  at  four  time  intervals:  a)  Day  45  (release);  b)  Day  60;  c)  Day  152;  d)  Day 
212  (end  of  simulation).  At  the  end  of  the  simulation,  the  highest  chlorophyll  concentrations 

are  in  mainstem. 
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8  Algae  as  Particles  -  Comparisons  with 
Observations 

Introduction 

Chlorophyll  concentrations  derived  from  the  model  of  algae  as  particles 
are  compared  to  observations  in  three  formats:  time  series,  seasonal 
spatial  distribution,  and  cumulative  distribution.  Concentrations  from  the 
conventional  model  with  algae  as  a  dissolved  substance  are  presented  in 
the  same  formats  for  comparison  with  the  alternate  model  and  with  the 
observations. 

Time  series  comparisons  are  presented  at  four  stations:  RET2.2,  RET2.4, 
LE2.2,  and  LE2.3  (Figure  8).  These  stations  encompass  the  spatial  extent 
of  the  spring  algal  bloom.  Chlorophyll  observations  were  collected  1  m 
below  the  surface  and  1  m  above  the  bottom  at  each  station,  once  or  twice 
monthly.  Samples  were  also  collected  1  m  above  and  1  m  below  the 
pycnocline  at  the  LE  stations.  These  samples  were  linearly  interpolated  to 
create  a  mid-depth  observation.  The  observations  are  compared  to  daily- 
average  model  outputs  at  station  locations  and  sample  depths.  The 
temporal  origin  of  the  plots  (Year  o)  is  the  commencement  of  the  model 
run  (November  1, 1993)  and  the  temporal  extent  is  one  year,  although  the 
run  extends  through  only  May  31. 

Seasonal  spatial  distributions  are  presented  for  two  seasons:  winter 
(December  -  February)  and  spring  (March  -  May).  The  mean  and  range  of 
surface  and  bottom  observations  were  computed  for  each  season  and 
compared  with  the  mean  and  range  of  the  daily  model  outputs.  The  spatial 
origin  of  the  plots  (km  o)  is  the  mouth  of  the  Potomac  Estuary  and  the 
spatial  extent  is  from  the  eastern  edge  of  the  computational  grid  (km  -10) 
to  the  Potomac  River  fall  line  (km  180).  The  preponderance  of  the  diatom 
bloom  occurs  below  km  120. 

The  cumulative  distribution  plots  summarize  and  compare  observations  at 
the  time  series  stations  with  model  results  at  the  same  locations  and  depths. 
The  observations  and  computations  are  sorted  from  smallest  to  largest 
values.  The  sorted  arrays  are  divided  into  quantiles  and  plotted  as  cumula¬ 
tive  distributions.  A  point  on  the  line  in  x-y  space  indicates  the  percentage 
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of  observations  or  computations  (x-axis)  less  than  the  indicated  concentra¬ 
tion  (y-axis).  The  50th  percentile  indicates  the  median  value.  Perfect 
correspondence  in  the  range  of  computed  and  observed  variables  is 
indicated  when  the  cumulative  distribution  of  modeled  values  exactly 
overlays  the  cumulative  distribution  of  observed  values. 

Time  Series  Comparisons 

The  temporal  characteristics  of  results  from  the  two  models  and  of 
comparisons  with  the  observations  vary  with  station  locations.  Some 
salient  characteristics  and  trends  can  be  discerned,  however.  The  particle 
model  demonstrates  high-frequency  oscillations,  with  periods  on  the  order 
of  days,  at  all  stations  (Figures  55  -  58).  The  amplitude  of  the  oscillations 
is  large  relative  to  the  magnitude  of  the  computations.  The  conventional 
model  demonstrates  comparable  oscillations  at  the  most  downstream 
station,  LE2.3  (Figure  55).  The  high-frequency  oscillations  are  damped 
with  distance  upstream,  however,  and  variations  with  longer  periods, 
weeks  to  months,  predominate  (Figures  56  -  58).  The  lengthy  intervals 
between  the  observations  confound  verification  of  the  oscillatory  behavior. 

At  the  most  downstream  station,  the  magnitude  of  results  from  both  models 
compares  with  the  observations  at  the  surface  and  mid-depth  while  both 
models  underestimate  the  chlorophyll  concentration  at  the  bottom 
(Figure  55).  Further  upstream,  the  particle  model  computes  lower  surface 
chlorophyll  concentrations  than  the  conventional  model  (Figures  56  -  58). 
The  magnitude  of  surface  chlorophyll  computed  by  the  particle  model 
reflects  the  magnitude  of  the  observations  better  than  the  conventional 
model.  Both  models  underestimate  observed  chlorophyll  at  the  bottoms  of 
the  LE  stations  (Figures  55, 56),  although  correct  magnitudes  are  computed 
by  both  models  at  the  RET  stations  (Figures  57,  58). 

Seasonal  Spatial  Distributions 

The  seasonal  spatial  distributions  “average  out”  the  high-frequency 
oscillations  in  the  computations  and  produce  more  readily  interpreted 
comparisons  between  the  characteristic  values  provided  by  the  models  and 
the  observations.  During  the  winter  months,  the  characteristics  of  the  two 
models  differ  substantially,  although  neither  is  apparently  superior  with 
regard  to  the  observations.  The  particle  model  computes  highest  surface 
concentrations  at  the  mouth  of  the  estuary;  concentrations  generally 
decrease  in  the  upstream  direction  and  are  lower  than  the  observations 
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Figure  55.  Observed  and  computed  chlorophyll  concentrations  at  station  LE2.3.  Results  are 
shown  from  the  model  with  algae  as  particles  at  (a)  surface,  (b)  mid-depth,  and  (c)  bottom 
and  from  a  conventional  model  at  (d)  surface,  (e)  mid-depth,  and  (f)  bottom. 
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Figure  56.  Observed  and  computed  chlorophyll  concentrations  at  station  LE2.2.  Results  are 
shown  from  the  model  with  algae  as  particles  at  (a)  surface,  (b)  mid-depth,  and  (c)  bottom 
and  from  a  conventional  model  at  (d)  surface,  (e)  mid-depth,  and  (f)  bottom. 
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Figure  57.  Observed  and  computed  chlorophyll  concentrations  at  station  RET2.4.  Results 
are  shown  from  the  model  with  algae  as  particles  at  (a)  surface  and  (b)  bottom  and  from  a 
conventional  model  at  (c)  surface  and  (d)  bottom. 
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Figure  58.  Observed  and  computed  chlorophyll  concentrations  at  station  RET2.2.  Results 
are  shown  from  the  model  with  algae  as  particles  at  (a)  surface  and  (b)  bottom  and  from  a 
conventional  model  at  (c)  surface  and  (d)  bottom. 
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between  km  30  and  80  (Figure  59).  The  conventional  model  computes 
highest  surface  computations  from  km  30  to  50  and  overestimates  the 
observations  in  this  reach  (Figure  59).  Peak  subsurface  computations  in 
both  models  occur  upstream  of  peak  surface  concentrations  and  are 
substantially  less  than  observations  in  both  cases  (Figure  59).  The  particle 
model  does  compute  a  subsurface  chlorophyll  maximum,  however,  circa  km 
40  to  60.  Mean  computed  bottom  concentrations  exceed  10  mg  nr3  while 
mean  computed  surface  concentrations  in  this  reach  never  exceed  6  mg  m-3. 


Figure  59.  Observed  and  computed  chlorophyll  concentrations  for  the  winter  season 
(December  -  February).  Km  zero  is  at  the  mouth  of  the  estuary.  Results  are  shown  from  the 
model  with  algae  as  particles  at  (a)  surface  and  (b)  bottom  and  from  a  conventional  model  at 

(c)  surface  and  (d)  bottom. 

The  shapes  of  the  spatial  distributions  from  both  models  are  preserved  from 
winter  to  spring.  The  particle  model  computes  maximum  surface  concentra¬ 
tions  at  the  mouth  of  the  estuary,  while  the  conventional  model  computes 
maximum  concentrations  between  km  30  and  50  (Figure  60).  The  particle 
model  provides  excellent  agreement  with  the  observations,  however,  while 
the  conventional  model  overestimates  surface  concentrations  at  all  stations 
(Figure  60).  Both  models  persist  in  locating  the  greatest  subsurface 
chlorophyll  concentrations  upstream  of  the  greatest  surface  concentrations. 
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Only  the  particle  model,  however,  computes  a  true  subsurface  maximum  in 
which  the  mean  chlorophyll  concentration  at  the  bottom,  ~  18  mg  nr3  circa 
km  50,  exceeds  the  mean  surface  concentration,  «  8  mg  m-3,  at  the  same 
location  (Figure  60). 


Potomac  River  (Part  Run37  BaseCase) 


Potomac  River  (R37  ICIWI  dupl) 


Figure  60.  Observed  and  computed  chlorophyll  concentrations  for  the  spring  season  (March  - 
May).  Km  zero  is  at  the  mouth  of  the  estuary.  Results  are  shown  from  the  model  with  algae  as 
particles  (a)  at  surface  and  (b)  bottom  and  from  a  (c)  conventional  model  at  surface  and  (d) 

bottom. 


Cumulative  Distribution  Plots 

Cumulative  distribution  plots  were  prepared  for  the  surface  observations 
and  corresponding  model  computations  (Figure  61),  for  the  bottom 
observations  and  corresponding  model  computations  (Figure  62),  and  for 
all  observations,  including  mid-depth,  and  corresponding  model  computa¬ 
tions  (Figure  63).  The  plot  at  the  surface  indicates  the  distribution  of  algae 
computed  as  particles  more  closely  represents  the  observed  distribution. 
The  preponderance  of  the  particle  distribution  is  below  the  observed,  while 
the  distribution  of  algae  as  dissolved  substance  exceeds  the  observations 
throughout.  Both  model  approaches  underestimate  the  distribution  of 
chlorophyll  observed  at  the  bottom.  A  noticeable  distinction  between  the 
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Figure  61.  Cumulative  distribution  plot  of  observed  and  modeled  surface  chlorophyll 

concentration. 


Figure  62.  Cumulative  distribution  plot  of  observed  and  modeled  bottom  chlorophyll 

concentration. 


ERDC/EL  TR-13-13 


64 


Percent  Less  Than 

Figure  63.  Cumulative  distribution  plot  of  observed  and  modeled  chlorophyll  concentration  at 

all  depths. 

two  approaches  is  that  half  the  particle  computations  are  exactly  zero,  while 
less  than  5%  of  the  dissolved  computations  are  zero.  The  cumulative 
distribution  plot,  which  considers  observations  and  comparable  computa¬ 
tions  at  all  stations  and  depths,  indicates  two  overarching  conclusions.  First, 
chlorophyll  computations  from  the  conventional  model  exceed  chlorophyll 
computations  from  the  particle  model.  The  median  concentration  from  the 
conventional  model  is  4.3  mg  m-3  versus  1.7  mg  m-3  from  the  particle  model. 
Second,  the  preponderance  of  observations  exceeds  both  models.  The 
median  observed  concentration  is  6.6  mg  nr3.  Both  models  compute 
extreme  concentrations,  which  exceed  the  maximum  observation,  however. 
This  effect  is  especially  noticeable  in  the  conventional  model  for  which  the 
upper  20%  of  the  chlorophyll  computations  exceeds  the  upper  20%  of  the 
observations. 

Time  Series  Analysis 

The  results  from  the  two  models  demonstrate  periodic  behavior.  The 
amplitudes  and  dominant  periods  vary  between  the  two  models,  as 
described  earlier  in  this  chapter.  The  periodic  behaviors  cannot  be 
validated  with  the  conventional  monitoring  data  due  to  the  lengthy  and 
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irregular  gaps  between  observations.  More  recent  data  from  alternate 
monitoring  programs  provides  insight  into  the  periodic  behavior  of  in-situ 
chlorophyll  and  allows  comparisons  with  the  model  patterns. 

Continuous  Monitoring  Data 

The  Chesapeake  Bay  Program  and  its  partner  agencies  commenced  a 
Shallow  Water  Monitoring  Program  circa  2005.  The  program  incorporates 
continuous  remote  sensors  in  locations  throughout  the  bay  system.  One 
station  was  operated  at  Nomini  Bay,  on  the  south  shore  of  the  Potomac 
adjacent  to  monitoring  station  LE2.2  (Figure  64),  from  2007  through 
2009.  Surface  chlorophyll  measures  (1.5  m  below  mean  low  water)  were 
recorded  at  15-minute  intervals  from  April  through  October  of  each  year. 
Observations  for  time  series  analysis  were  retrieved  from  an  on-line  data¬ 
base  (VIMS  2012)  for  the  period  April  -  June  2008.  Although  the 
observations  were  not  concurrent  with  the  model  application,  analyses 
were  conducted  based  on  the  assumption  that  characteristic  periodic 
behaviors  are  consistent  from  year  to  year. 

A  second  station  was  operated  at  St.  George’s  Island,  on  the  northern  shore 
of  the  Potomac  (Figure  64),  from  2006  through  2012.  Surface  chlorophyll 
measures  (1.3  m  below  mean  low  water)  were  recorded  at  10  to  15  minute 
intervals  from  April  through  October  of  these  years.  Inspection  of  the  data 
(MD  DNR  2012)  indicated  the  record  from  mid-April  through  June  2012 
provided  the  best  opportunity  for  time  series  analysis  due  to  regular 
sampling  intervals  and  limited  data  gaps.  Model  results  for  comparison  with 
these  observations  were  retrieved  from  the  surface  cell  at  Station  LE2.2  for 
the  period  January  1  through  May  31, 1994. 

The  observations  (Figure  65,  66)  demonstrate  three  properties:  regular 
oscillations  with  a  period  on  the  order  of  days,  longer-term  variations  with 
periods  of  tens  of  days,  and  irregular  large-amplitude  spikes.  Qualitatively, 
the  observations  resemble  the  results  from  the  model  with  algae  as  a 
dissolved  substance  (Figure  67),  except  the  model  lacks  the  spikes.  The 
spikes  and  regular  oscillations  are  present  in  the  particle  model  (Figure  68), 
but  the  longer-term  variations  are  difficult  to  perceive.  The  particle  model  is 
distinctive  in  that  it  computes  concentrations  of  zero,  which  are  not  present 
in  either  the  observations  or  the  dissolved  model. 
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Figure  64.  Location  of  Nomini  Bay  and  St.  George’s  Island  continuous  monitoring  stations. 
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Figure  65.  Observed  chlorophyll  at  Nomini  Bay,  April  -  May  2008. 


Figure  66.  Observed  chlorophyll  at  St.  George’s  Island,  April  -  May  2012. 
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Figure  67.  Surface  chlorophyll  computed  at  Station  LE2.2  using  a  conventional  model  of 

algae  as  a  dissolved  substance. 


Figure  68.  Surface  chlorophyll  computed  at  Station  LE2.2  using  a  model  of  algae  as  particles. 
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Harmonic  Analysis 

A  linear  model  was  proposed  to  describe  the  observed  record  and  the 
predicted  time  series  at  LE2.2: 


chl  =  a  +  bl*t  +  b2*sin 


2nt 

12.4 


b3*cos 


2nt 

12.4 


+  b4*sin 


2nt 

24 


b5*cos 


2nt 

24 


(9) 


in  which: 

chi  =  Chlorophyll  concentration  (mg  m  3) 
t  =  Time 

a,  bi ...  bs  =  Model  parameters 

The  b2  and  b3  parameters  represent  the  amplitudes  of  harmonics  with  the 
period  of  the  lunar  semi-diurnal  tide,  12.4  hours.  The  b4  and  bs  parameters 
represent  harmonics  with  a  daily  period.  The  amplitudes  of  the  sin  and  cos 
can  be  combined  using  the  identity: 

A  sin  (a)  +  B  cos  (a)  =  >/a2  +B2  sin|a  +  tan_1^^Jj  (10) 

Parameters  for  the  observations  and  the  predicted  time  series  were  deter¬ 
mined  via  linear  regression.  The  fraction  of  variation  explained  by  the 
models  was  typically  small  (R2  ~  0.1)  but  highly  significant  (p  <  0.001).  The 
amplitudes  of  the  semi-diurnal  and  diurnal  harmonics  were  roughly  equal 
for  the  observations  and  for  the  model  of  algae  as  dissolved  substance 
(Figure  69).  The  influences  of  tides  and  of  diurnal  irradiance  were  equiva¬ 
lent.  The  semi-diurnal  amplitude  for  the  particle  model  was  several  times 
larger  than  any  other  amplitude,  while  the  diurnal  amplitude  was  the  least 
of  all.  The  particle  model  exhibited  a  strong  influence  from  tides  but  almost 
no  influence  from  diurnal  irradiance. 

Power  Spectra 

Power  spectra  for  the  observed  and  computed  time  series  were  obtained 
using  the  Fast  Fourier  Transform  feature  of  MATLAB  (MathWorks  2012). 
Power  spectra  of  the  observed  time  series  indicate  substantial  peaks  at  the 
semi-diurnal  and  diurnal  periods.  The  Nomini  Bay  spectrum  (Figure  70) 
indicates  more  power  in  the  semi-diurnal  period.  The  St.  Georges’s 
spectrum  (Figure  71)  indicates  more  power  in  the  diurnal  period.  The 
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Nomini  St  Georges  LE2*2  particle  LE2.2  dissolved 


Figure  69.  Amplitudes  of  lunar  semi-diurnal  and  diurnal  harmonics  of  observed  and 
computed  chlorophyll  time  series.  Observations  from  Nomini  Bay  and  St.  George's  island. 
Model  results  at  LE2.2  using  algae  as  particles  and  as  dissolved  substance. 


Period  (days) 


Figure  70.  Power  spectrum  of  chlorophyll  observations  at  Nomini  Bay.  The  red  arrows  indicate 

the  semi-diurnal  and  diurnal  periods. 
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Figure  71.  Power  spectrum  of  chlorophyll  observations  at  St.  George's  Island.  The  red  arrows 
indicate  the  semi-diurnal  and  diurnal  periods. 


relative  magnitude  of  these  two  peaks  corresponds  to  the  relative  magnitude 
of  the  amplitudes  in  the  harmonic  analysis  (Figure  69).  Both  observed 
spectra  exhibit  power  roughly  equivalent  to  the  semi-diurnal  and  diurnal 
periods  at  periods  of  seven  days  or  more. 

The  spectrum  from  the  dissolved  model  (Figure  72)  resembles  the 
observations  by  demonstrating  peaks  with  nearly  equivalent  power  at  the 
semi-diurnal  and  diurnal  periods.  As  with  the  observations,  power  nearly 
equivalent  to  the  semi-diurnal  and  diurnal  periods  occurs  at  periods  greater 
than  10  days.  The  particle  model  (Figure  73)  demonstrates  its  greatest 
power  at  the  semi-diurnal  period.  Diurnal  fluctuations  are  virtually  absent. 
The  power  in  both  models  at  periods  greater  than  10  days  is  equivalent.  The 
particle  model  has  orders  of  magnitude  more  power  than  the  dissolved 
model  at  high  frequencies,  less  than  0.5  days.  These  high-frequency 
fluctuations  are  also  present  at  Nomini  Bay  but  not  St.  George’s. 
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Period  (days) 


Figure  72.  Power  spectrum  of  surface  chlorophyll  computed  at  Station  LE2.2  using  a 
conventional  model  of  algae  as  a  dissolved  substance.  The  red  arrows  indicate  the  semi¬ 
diurnal  and  diurnal  periods.  Note  the  equivalent  power  at  these  two  periods. 


Figure  73.  Power  spectrum  of  surface  chlorophyll  computed  at  Station  LE2.2  using  the  model 
of  algae  as  particles.  The  red  arrows  indicate  the  semi-diurnal  and  diurnal  periods.  Note  the 
predominance  of  the  semi-diurnal  period  and  the  absence  of  power  at  the  diurnal  period. 
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Summary 

Results  from  the  model  with  algae  treated  as  particles  were  compared  to 

observations  and  to  results  from  a  conventional  model,  with  algae  as 

dissolved  substance,  in  a  variety  of  formats.  Salient  results  are: 

•  The  model  with  algae  as  particles  is  superior  to  the  conventional  model 
in  two  regards.  The  particle  model  computes  high  chlorophyll 
concentrations  at  great  depths  and  the  particle  model  produces  a  true 
subsurface  chlorophyll  maximum. 

•  The  spatial  distributions  of  chlorophyll  computed  by  the  two  models 
differ.  The  particle  model  computes  maximum  chlorophyll  at  the 
mouth  of  the  estuary.  The  conventional  model  computes  maximum 
chlorophyll  further  upstream.  Surface  chlorophyll  concentrations 
computed  by  the  conventional  model  are  greater  than  the  particle 
model.  Neither  model  is  clearly  superior  with  regard  to  reproducing  the 
observations. 

•  Both  models  underestimate  the  cumulative  distribution  of  the 
observations.  The  distribution  from  the  conventional  model  is  closer  to 
the  distribution  from  the  observations,  since  the  excess  chlorophyll 
computed  at  the  surface  partially  compensates  for  the  absence  of 
computed  chlorophyll  in  deep  water  near  the  mouth  of  the  estuary. 

•  The  particle  model  exhibits  high-frequency  (period  <  l  day)  oscillations 
and  spikes  that  are  absent  from  the  model  with  algae  as  dissolved 
substance.  The  oscillations  and  spikes  are  present  in  observations  from 
at  least  one  station  with  a  database  sufficient  for  analysis.  Power 
spectra  of  observations  and  of  computations  from  the  dissolved  model 
indicate  significant  peaks  at  the  semi-diurnal  (tidal)  and  diurnal  (daily) 
periods.  The  power  spectrum  of  the  particle  model  shows  its  greatest 
power  at  the  semi-diurnal  period,  indicating  strong  influence  from 
tidal  transport.  Diurnal  fluctuations,  produced  by  the  daily  cycle  of 
solar  irradiance,  are  completely  absent,  however. 

•  Results  from  the  particle  model  are  “patchy”  in  time  and  space.  Nearly 
40%  of  the  computations  corresponding  to  observations  are  zero.  The 
large  fraction  of  null  concentrations  does  not  correspond  to 
observations. 

•  Results  from  the  conventional  model  reflect  intense  effort  to  tune  the 
model  to  observations.  For  comparison  purposes,  the  particle  model 
used  an  identical  parameter  set.  Improved  results  may  be  attainable 
with  the  particle  model  following  development  of  a  specific  parameter 
set. 
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9  Algae  as  Particles  -  Sensitivity  Analyses 

The  basic  characteristics  of  particle  transport  and  of  algae  modeled  as 
particles  have  been  established  in  preceding  chapters.  The  present  chapter 
examines  sensitivity  of  particles  and  algae  to  basic  particle  behaviors  and 
then  investigates  the  potential  effects  of  different  hydrology  on  the 
previous  conclusions. 

Sensitivity  to  Particle  Behavior 

Chapter  1  describes  a  conceptual  model  of  the  Potomac  River  spring 
diatom  bloom.  The  model  proposes  that  algae  regulate  their  buoyancy  to 
maintain  their  residence  in  the  central  portion  of  the  estuary.  The  authors 
perform  three  experiments  here  that  initiate  the  investigation  into  the 
potential  effects  of  buoyancy  regulation  on  particle  residence  time  and 
location  and  on  chlorophyll  concentration: 

1.  Particles  rise  at  the  rate  of  l  m  d1. 

2.  Particles  sink  at  the  rate  of  l  m  d-1. 

3.  Particles  sink  at  the  rate  of  1  m  d1  in  the  LE  portion  of  the  estuary  and  rise 
at  the  rate  1  m  d1  in  the  RET  portion  of  the  estuary. 

For  the  examination  of  particle  residence  time  and  location,  the  kinetics 
rules  that  keep  the  number  of  particles  in  the  system  constant  are 
disabled.  Particles  are  free  to  exit  the  system  or  to  adhere  to  the  bottom. 
The  kinetics  rules  remain  in  place  for  the  examination  of  chlorophyll.  The 
three  sensitivity  runs  are  compared  to  the  base  particle  simulation  and  to 
the  base  run  of  algae  as  particles. 

Particle  residence  time  and  location 

The  tracker  statistics  provide  quantitative  measures  of  the  particle  residence 
time  and  location.  The  initial  distribution  of  particles  (Figure  74a  -  74d)  is 
the  same  in  all  runs  and  as  described  in  Chapter  6: 10,000  particles  are 
released  between  km  30  and  80  on  Day  45  of  the  model  run.  By  Day  60, 

15  days  after  the  particle  release,  the  run  with  buoyant  particles  indicates 
that  particles  move  into  the  surface  layer  and  downstream  relative  to  the 
base  condition  with  neutrally-buoyant  particles  (Figure  75a,  75b).  The  run 
in  which  particles  sink  (Figure  75c)  and  the  run  in  which  particles  rise  and 
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Figure  74.  Initial  particle  distribution  (Day  45)  for  four  model  runs:  a)  base  case;  b)  particles 
rise  at  1  m  d* 1;  c)  particles  sink  at  1  m  d1;  and  d)  particles  rise  in  RET  and  sink  in  LE. 


Figure  75.  Particle  distribution  at  Day  60  for  four  model  runs:  a)  base  case;  b)  particles  rise  at 

1  m  d1;  c)  particles  sink  at  1  m  d1;  and  d)  particles  rise  in  RET  and  sink  in  LE. 
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sink  (Figure  75d)  are  similar  to  each  other.  Particle  concentration 
diminishes  in  the  surface  layer  of  the  LE  region  and  increases  in  the 
subsurface  layers  of  the  RET  reach.  The  presence  of  particles  in  the 
subsurface  layers  of  RET  while  particles  are  rising  indicates  transport 
upstream  from  the  lower  estuary,  in  keeping  with  the  conceptual  model. 

By  Day  150,  most  of  the  buoyant  particles  have  left  the  system  (Figure  76b). 
The  remaining  fraction  is  substantially  lower  than  the  remaining  fraction  of 
neutrally-buoyant  particles,  indicating  a  characteristic  of  particle  buoyancy. 
Particles  which  rise  tend  to  leave  the  system  on  the  net  outflowing  surface 
currents.  Few  particles  have  left  the  system  in  the  runs  in  which  particles 
sink  (Figure  76c)  and  in  which  particles  rise  and  sink  (Figure  76d).  Instead, 
particles  accumulate  in  the  mid-depth  and  bottom  of  the  RET  reach.  These 
runs  indicate  that  sinking  tends  to  increase  the  residence  time  of  particles  in 
the  system. 


Figure  76.  Particle  distribution  at  Day  150  for  four  model  runs:  a)  base  case;  b)  particles  rise 
at  1  m  d1;  c)  particles  sink  at  1  m  d1;  and  d)  particles  rise  in  RET  and  sink  in  LE. 


The  characteristics  noted  at  Day  150  persist  at  the  completion  of  the  runs 
(Figure  77).  The  run  with  buoyant  particles  retains  fewer  particles  than  the 
run  with  neutrally  buoyant  particles.  The  runs  in  which  particles  sink 
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accumulate  particles  in  the  subsurface  layers  of  the  RET  reach.  The  run  in 
which  particles  rise  and  sink  (Figure  76d)  is  not  materially  different  from 
the  run  in  which  particles  sink  exclusively  (Figure  76c).  Apparently,  the 
buoyancy  is  not  sufficient  to  raise  particles  to  the  surface.  The  particles  are 
effectively  trapped  at  the  bottom. 


Figure  77.  Particle  distribution  at  run  completion  (Day  210)  for  four  model  runs:  a)  base  case;  b) 
particles  rise  at  1  m  d1;  c)  particles  sink  at  1  m  d1;  and  d)  particles  rise  in  RET  and  sink  in  LE. 


Chlorophyll  Concentration 

Chlorophyll  concentrations  are  examined  using  the  seasonal-average 
longitudinal  plots.  As  expected,  particle  buoyancy  diminishes  the  chloro¬ 
phyll  concentration  along  the  bottom  in  both  winter  (Figure  78b)  and  spring 
(Figure  80b)  relative  to  the  run  with  neutrally  buoyant  particles.  Particle 
sinking  increases  the  amount  of  chlorophyll  at  the  bottom  upstream  of  km 
50  (Figure78  c  and  d,  Figure  80  c  and  d).  The  run  in  which  particles  sink 
exclusively  increases  bottom  chlorophyll  concentration  upstream  as  far  as 
km  100.  The  run  with  rising  and  sinking  demonstrates  the  influence  of 
buoyancy  at  the  head  of  the  algal  bloom,  circa  km  80  to  100.  Algae  in  this 
reach  apparently  rise  to  the  surface.  No  run  matches  the  elevated  observed 
chlorophyll  concentrations  in  the  lower  40  km  of  the  estuary. 
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Figure  78.  Observed  and  computed  winter,  bottom,  chlorophyll  concentrations  for  four  model 
runs:  a)  base  case;  b)  particles  rise  at  1  m  d1;  c)  particles  sink  at  1  m  d1;  and  d)  particles  rise 

in  RET  and  sink  in  LE. 


Potomac  River  (Part  Run39) 
Surface  Chlorophyll  Winter  1994 


Potomac  River  (Part  Run40) 
Surface  Chlorophyll  Winter  1994 


Figure  79.  Observed  and  computed  winter,  surface,  chlorophyll  concentrations  for  four  model 
runs:  a)  base  case;  b)  particles  rise  at  1  m  d1;  c)  particles  sink  at  1  m  d1;  and  d)  particles  rise 

in  RET  and  sink  in  LE. 
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Figure  80.  Observed  and  computed  spring,  bottom,  chlorophyll  concentrations  for  four  model 
runs:  a)  base  case;  b)  particles  rise  at  1  m  d1;  c)  particles  sink  at  1  m  d1;  and  d)  particles  rise 

in  RET  and  sink  in  LE. 


The  run  with  buoyant  particles  provides  superior  surface  chlorophyll 
computations  in  winter  (Figure  79b)  and  computations  equivalent  to  the 
neutrally  buoyant  run  in  spring  (Figure  81b).  At  the  surface,  the  run  in 
which  particles  sink  exclusively  provides  the  least  satisfactory  observations 
(Figures  79c,  81c).  The  run  in  which  particles  rise  and  sink  provides  a 
surface  chlorophyll  peak  circa  km  90  and  partially  validates  the  conceptual 
model,  especially  during  winter  (Figure  79d),  although  elsewhere  this 
model  shares  the  unsatisfactory  properties  of  the  run  in  which  particles 
sink  exclusively. 

Chlorophyll  Distribution  Plot 

The  distributions  of  chlorophyll  computed  by  the  base  model,  the  model  in 
which  particles  sink,  and  the  model  in  which  particles  rise  and  sink  are 
virtually  congruent  through  75%  of  their  distributions  (Figure  82).  The 
chlorophyll  distribution  in  the  run  with  buoyant  particles  is  lower  than  the 
others  through  90%  of  its  distribution.  The  differences  in  the  distributions 
are  primarily  in  the  extreme  values.  The  run  in  which  particles  sink 
exclusively  provides  the  larger  portion  of  extreme  values,  indicating  the 
effects  of  particle  accumulation  in  a  few  patches.  The  chlorophyll  distribu¬ 
tions  from  all  runs  are  characteristically  less  than  the  observations. 
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Figure  81.  Observed  and  computed  spring,  surface,  chlorophyll  concentrations  for  four  model 
runs:  a)  base  case;  b)  particles  rise  at  1  m  d1;  c)  particles  sink  at  1  m  d1;  and  d)  particles  rise 

in  RET  and  sink  in  LE. 


Figure  82.  Cumulative  distribution  plot  of  observed  and  computed  chlorophyll  concentration 
at  all  depths  and  stations.  Computed  results  from  four  model  runs:  base  case;  particles  rise 
at  1  m  d1;  particles  sink  at  1  m  d -1;  and  particles  rise  in  RET  and  sink  in  LE. 
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Sensitivity  to  Hydrology 

The  year  1994  selected  for  the  base  runs  is  a  year  of  average  hydrology  in  the 
Chesapeake  Bay  system  and  of  above-average  runoff  in  the  Potomac  River. 
The  year  1999  is,  by  contrast,  a  dry  year  in  both  the  Chesapeake  Bay  system 
and  in  the  Potomac.  Average  flow  at  the  Potomac  fall  line  for  the  period 
November  1998  through  May  1999  is  208  m3  s-1,  compared  to  802  m3  s-1  for 
the  period  November  1993  through  May  1994.  The  1999  period  lacks  the 
number  of  flood  events  and  the  peak  flows  present  in  1994  (Figure  83).  To 
test  the  sensitivity  of  model  results  to  hydrology,  the  model  of  algae  as 
particles  was  run  for  the  period  November  1998  through  May  1999,  using 
the  same  initial  conditions  and  parameter  set  as  the  base  run.  The  tracker 
statistics  indicate  the  effects  of  the  dry  hydrology  as  early  as  15  days  after 
the  particle  release  (Figure  84b,  f).  More  particles  penetrate  further 
upstream  in  1999  than  in  1994.  This  characteristic  persists  through  the 
remainder  of  the  run.  The  dry  year  shows  more  particles  upstream  of  km  60 
than  the  wet  year.  Since  the  number  of  particles  is  constant,  the  wet  year  has 
more  particles  in  the  downstream  portion  of  the  estuary  than  the  dry  year. 
The  change  in  spatial  distribution  of  particles  is  entirely  consistent  with  the 
relative  magnitudes  of  the  net  flows  in  the  two  years. 


Figure  83.  Fall-line  flow  in  the  Potomac  River  for  dry  (1999)  versus  wet  (1994)  hydrologic 

conditions. 
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Figure  84.  Particle  distribution  computed  for  dry  hydrology  (1999)  versus  base  (1994) 
conditions.  Base  conditions  are  in  the  left  column  (panels  a  -  d).  The  dry  conditions  are  in  the 

right  column  (panels  e  -  h). 
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The  periodic  behavior  and  the  vertical  distributions  of  computed  chloro¬ 
phyll  are  similar  in  the  two  years  (Figure  85).  Both  runs  demonstrate  high- 
frequency  oscillations  and  both  runs  reflect  the  magnitude  of  observed 
surface  chlorophyll  while  underestimating  the  observed  chlorophyll  in  the 
deep  waters  of  the  lower  estuary.  A  large  algal  bloom  occurred  in  April  1999 
(0.5  years  from  model  commencement)  that  is  not  represented  in  the 
simulation. 
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Figure  85.  Observed  and  computed  chlorophyll  at  LE2.2  for  three  depths  (surface,  mid-depth, 
bottom)  and  two  hydrologic  conditions:  1994,  wet  (panels  a  -  c)  and  1999,  dry  (panels  d  -f). 
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Despite  the  difference  in  hydrology,  the  computed  spatial  chlorophyll 
distributions  are  similar  in  the  two  years  (Figure  86).  Computed  chlorophyll 
in  the  bottom  is  restricted  to  the  reach  between  km  30  and  km  110.  No 
bottom  chlorophyll  is  computed  in  the  lower  30  km  of  the  estuary  despite 
the  high  concentrations  observed.  The  computed  surface  chlorophyll  is 
shifted  farther  downstream  in  spring  1994  compared  to  1999,  perhaps 
driven  by  the  higher  runoff.  Magnitudes  of  computed  surface  chlorophyll 
are  similar  in  the  two  years,  however.  Model  results  in  both  years  match  the 
surface  observations  with  the  exception  of  the  bloom  event  centered  at  km 
30  during  1999. 


Figure  86.  Observed  and  computed  spring  chlorophyll  concentrations  for  two  hydrologic 
conditions:  1994,  wet  (panels  a,  b)  and  1999,  dry  (panels  c,  d). 


Comparisons  of  observed  and  computed  chlorophyll  distributions  in  1999 
(Figure  87)  resemble  the  comparisons  from  1994  (Figure  63).  The  com¬ 
puted  distribution  is  less  than  observed  and  demonstrates  a  substantial 
fraction  of  zero  values  which  are  not  present  in  the  observations.  The 
chlorophyll  results  from  1999  indicate  the  model  characteristics  demon¬ 
strated  in  the  base  run  can  be  generalized  across  varying  hydrologic 
conditions.  The  large  alteration  in  hydrology  from  1994  to  1999  induced  no 
aberrant  behavior  in  the  model  and  the  model  is  suitable  for  further 
investigation  and  application. 
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Figure  87.  Cumulative  distribution  plot  of  observed  and  computed  chlorophyll  concentration 

at  all  depths  and  stations  for  1999. 
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10  Summary  and  Conclusions 

Introduction 

Reliable  prediction  of  phytoplankton  transport  and  production  is  central  to 
the  understanding  and  remediation  of  a  host  of  environmental  problems, 
including  eutrophication  and  the  occurrence  of  harmful  algal  blooms. 
Quantitative  phytoplankton  models  have  existed  for  decades,  yet  virtually 
all  examples  are  based  on  the  original  framework  wherein  phytoplankton 
are  a  dissolved  substance  transported  passively  and  exclusively  by 
hydrodynamic  processes.  Herein,  the  authors  test  the  hypothesis  that 
phytoplankton  dynamics  (particularly  the  occurrence  of  blooms)  can  be 
more  accurately  predicted  by  treating  phytoplankton  as  discrete  particles 
capable  of  self-induced  transport  via  buoyancy  regulation  or  other 
behaviors. 

ICM  Particle  Tracking 

The  transport  of  discrete  particles  is  accommodated  through  development 
of  a  particle-tracking  algorithm  which  works  within  the  CE-QUAL-ICM 
eutrophication  model.  Incorporation  of  the  algorithm  into  ICM  is  necessary 
in  order  for  particles  to  respond  to  environmental  cues  computed  by  the 
eutrophication  model.  The  calculations  are  performed  within  the  physical 
plane  through  the  application  of  general  algebra,  geometry,  and  physics. 

The  effects  of  turbulent  diffusion  and  dispersion  are  translated  into  random 
displacement  in  the  particle  model.  The  random  displacement  is  added  to 
the  deterministic  location  derived  from  advection.  As  a  first  approach, 
spatially  uniform  diffusion  and  dispersion  are  considered. 

Phytoplankton  Kinetics 

Algae  are  quantified  as  biomass  (carbon)  per  particle.  The  number  of 
particles  in  the  system  is  constant.  Changes  in  algal  biomass  are  represented 
by  altering  the  biomass  attached  to  a  particle  rather  than  by  changing  the 
number  of  particles.  The  particle-based  kinetics  are  derived  from  phyto¬ 
plankton  kinetics  in  the  CE-QUAL-ICM  eutrophication  model. 

The  rate  of  change  in  biomass  attached  to  a  particle  is  determined  by  the 
environment  surrounding  the  particle;  e.g.,  light  and  temperature.  The  rate 
of  change  also  influences  the  environment;  e.g.,  ambient  nutrient 
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concentration.  A  means  is  necessary  to  translate  from  the  particle-based 
system  to  the  continuum.  The  approach  employed  here  takes  advantage  of 
the  discrete  computational  grid  employed  by  the  ICM  eutrophication  model. 
The  particles  are  arrayed  in  a  computational  structure  that  labels  each 
particle  with  its  current  grid  cell.  The  structure  and  labeling  allow  the 
computer  code  to  loop  over  all  particles  or  over  all  cells.  The  code  can  query 
particles  “What  cell  am  I  in?”  or  it  can  query  cells  “What  particles  do  I 
contain?”  The  revised  ICM  code  with  particle-based  kinetics  employs  the 
algal  subroutine  from  the  conventional  ICM  code  to  provide  the  particle 
environment,  compute  growth  and  loss  rates,  and  to  provide  interactions 
with  quantities  computed  in  the  continuum  such  as  nutrient  concentration. 

Particle-tracking  models  conventionally  incorporate  rules  for  particle 
behavior.  Rules  are  especially  common  for  determining  particle  behavior 
when  the  trajectory  takes  a  particle  outside  the  model  domain.  Three  rules 
are  enforced  that  keep  the  number  of  particles  in  the  system  constant  and 
ensure  that  each  particle  is  associated  with  viable  algal  biomass.  Particles 
that  exit  through  an  open  boundary,  that  adhere  to  the  bottom,  or  that 
contain  no  viable  algae  are  eliminated  and  replaced  by  splitting  the 
particle  containing  the  largest  algal  biomass  in  two. 

The  Potomac  River  Estuary 

The  Potomac  River  Estuary  is  a  major  sub-estuary  of  the  larger 
Chesapeake  Bay  system.  The  estuary  extends  190  km  from  the  junction 
with  Chesapeake  Bay  to  the  head  of  tide  at  Washington  DC.  Mean  tide 
range  near  the  mouth  is  0.38  m  and  is  0.84  m  at  Washington  DC.  The 
estuary  is  a  drowned  river  valley  and  the  saline  portion  is  weakly  stratified. 
The  primary  freshwater  source  is  from  the  29,940  km2  upland  watershed 
and  enters  at  the  head  of  tide.  Lesser  volumes  enter  the  estuary  from  the 
adjacent  watershed  below  the  head  of  tide.  Long-term  mean  runoff  at 
Washington  DC  is  339  m3  s_1.  This  long-term  mean  is  subject  to  regular 
seasonal  fluctuations  as  well  as  extremes  due  to  flood  and  drought. 

The  particle-tracking  model  is  inserted  into  the  Chesapeake  Bay 
Environmental  Model  Package  (CBEMP).  The  CBEMP  consists  of  three 
independent  models:  a  watershed  model  (WSM),  a  hydrodynamic  model 
(HM),  and  a  eutrophication  model  (WQM).  The  WSM  provides  distributed 
flows  to  the  HM  and  nutrient  and  solids  loads  to  the  WQM.  The  HM 
computes  three-dimensional  intra-tidal  transport  and  supplies  transport 
parameters  to  the  WQM  on  an  hourly  basis.  The  WQM  computes  algal 
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biomass,  nutrient  cycling,  and  dissolved  oxygen,  as  well  as  numerous 
additional  constituents  and  processes. 

Both  the  HM  and  WQM  operate  on  a  three-dimensional  grid  which 
encompasses  the  entire  Chesapeake  Bay  system.  For  this  study,  the 
Potomac  River  portion  of  the  grid  was  extracted  along  with  the  associated 
hydrodynamics.  An  adjacent  portion  of  Chesapeake  Bay  was  included  so 
that  downstream  boundary  conditions  could  be  specified  at  a  distance 
sufficient  to  minimize  influence  on  the  upstream  portions  of  the  system. 
The  resulting  grid  extended  the  190  km  length  of  the  Potomac  River 
estuary  and  incorporated  numerous  embayments  and  tributaries.  The 
extracted  model  was  checked  against  the  original  model  and  against  data 
to  ensure  that  the  extraction  was  performed  correctly.  This  model  version 
was  used  as  a  starting  point  for  the  particle-tracking  algorithms  and 
associated  algal  kinetics. 

Basic  Particle  Simulation 

The  basic  particle  simulation  was  determined  following  a  variety  of 
numerical  experiments  and  sensitivity  runs.  The  base  simulation  extends 
from  November  1, 1993  through  May  31, 1994.  The  duration  is  selected  to 
encompass  the  spring  algal  bloom  period.  The  year  1994  corresponds  to  an 
average  hydrologic  year  selected  for  analysis  of  Chesapeake  Bay  model 
results,  although  flows  in  the  Potomac  during  the  interval  are  above- 
average.  Ten-thousand  particles  are  released  45  days  into  the  run  and  are 
spaced  uniformly  from  top-to-bottom  in  the  channel  of  the  river  between 
km  30  and  80.  Particles  are  distributed  in  model  grid  cells  with  the 
number  in  each  cell  specified  to  yield  4  x  io-6  particles  nr3.  The  release 
timing  and  region  correspond  to  the  typical  appearance  date  and  expanse 
of  bloom  algae. 

The  experiments  with  particle  tracking  indicate  that  particle  transport  is 
fundamentally  different  than  transport  of  a  dissolved  substance.  In 
particular,  the  residence  time  of  particles  released  in  the  center  of  the 
Potomac  Estuary  is  longer  than  the  residence  time  of  a  comparable  release 
of  dissolved  substance.  Ninety  days  after  the  release,  65%  of  the  particles 
remain  in  the  model  domain  while  the  dissolved  substance  has  disappeared. 
A  second  distinction  is  in  the  “patchy”  distribution  of  particles  compared  to 
the  continuous  distribution  of  dissolved  substance. 
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Algae  as  Particles  -  Spatial  and  Temporal  Distribution 

The  base  model  run  with  algae  as  particles  accompanies  the  basic  particle 
simulation.  Each  particle  is  initiated  with  2.5  x  ios  g  algal  carbon,  yielding 
1.0  g  C  m-3  at  initiation  of  the  simulation.  Parameters  for  the  algal  kinetics 
are  drawn  from  the  spring  diatom  group  of  the  CBEMP.  These  values  are 
employed  to  provide  comparison  to  algae  as  simulated  in  conventional 
models.  No  attempt  is  made,  at  this  stage  of  development,  to  optimize  the 
performance  of  the  model  with  algae  as  particles. 

A  conventional  model  run  was  completed  to  provide  comparison  with  the 
particle  approach.  Algae  were  initiated  at  a  concentration  of  1.0  g  C  m-3  at 
the  same  time  and  in  the  same  extent  of  the  Potomac  River  as  the  particle 
release.  Concentration  boundary  conditions  at  the  open  ends  of  the  model 
domain  were  specified  as  zero.  Kinetics  parameters  were  identical  to  the 
particle  simulation  with  one  exception.  The  conventional  run  employed  an 
algal  settling  velocity,  as  in  the  Chesapeake  Bay  model.  Settling  was 
omitted  from  the  base  particle  simulation. 

Particle  model  results  rendered  as  elevation  views  along  the  river  transect 
demonstrate  a  heterogeneous  spatial  distribution,  which  manifests 
immediately  upon  release  of  the  particles.  Chlorophyll  concentrations  near 
the  bottom  commonly  exceed  concentrations  at  lesser  depths  above.  The 
surface  layer,  which  receives  the  greatest  illumination,  is  at  times  devoid  of 
algae.  Chlorophyll  concentrations  computed  in  the  conventional  model  are 
much  more  smoothly  distributed.  Soon  after  bloom  initiation,  there  is  an 
indication  of  a  subsurface  chlorophyll  maximum  caused  by  upstream 
advection  of  the  initial  uniform  chlorophyll  distribution.  Otherwise,  the 
vertical  chlorophyll  distribution  declines  monotonically  from  the 
illuminated  surface  waters  to  the  bottom.  No  surface  cells  within  the 
bloom  region  are  devoid  of  algae. 

Surface  chlorophyll  concentrations  computed  by  the  particle  model  follow 
the  heterogeneous  pattern  observed  in  the  elevation  views.  Chlorophyll 
penetrates  into  the  smallest  tributaries  and,  by  the  end  of  the  simulation, 
the  highest  concentrations  are  in  the  tributaries  rather  than  in  the  surface 
waters  of  the  mainstem  Potomac.  The  conventional  model  provides  more 
uniform  surface  chlorophyll  distributions.  At  the  end  of  the  simulation, 
chlorophyll  concentrations  in  the  mainstem  are  as  high,  or  higher,  than  in 
the  tributaries.  The  contrast  between  the  uniform  surface  distribution  in 
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the  conventional  model  and  the  patchy  distribution  in  the  particle  model  is 
striking. 

Algae  as  Particles  -  Comparisons  with  Observations 

Results  from  the  model  with  algae  treated  as  particles  were  compared  to 
observations  and  to  results  from  a  conventional  model,  with  algae  as 
dissolved  substance,  in  a  variety  of  formats.  Salient  results  are  the  following: 

•  The  model  with  algae  as  particles  is  superior  to  the  conventional  model 
in  two  regards.  The  particle  model  computes  high  chlorophyll 
concentrations  at  great  depths  and  the  particle  model  produces  a  true 
subsurface  chlorophyll  maximum. 

•  The  spatial  distributions  of  chlorophyll  computed  by  the  two  models 
differ.  Neither  model  is  clearly  superior  with  regard  to  reproducing  the 
observations. 

•  The  particle  model  exhibits  high-frequency  (period  <  l  day)  oscillations 
and  spikes  that  are  absent  from  the  model  with  algae  as  dissolved 
substance.  The  oscillations  and  spikes  are  present  in  observations  from 
at  least  one  station  with  a  database  sufficient  for  analysis.  The  power 
spectrum  of  the  particle  model  shows  its  greatest  power  at  the  semi¬ 
diurnal  period,  indicating  strong  influence  from  tidal  transport. 

Diurnal  fluctuations,  produced  by  the  daily  cycle  of  solar  irradiance, 
are  completely  absent,  however. 

•  Results  from  the  particle  model  are  “patchy”  in  time  and  space.  Nearly 
40%  of  the  computations  corresponding  to  observations  are  zero.  The 
large  fraction  of  null  concentrations  does  not  correspond  to 
observations. 

Algae  as  Particles  -  Sensitivity  Analyses 

Two  sensitivity  analyses  were  performed.  The  first  examined  the 
sensitivity  of  particle  location  and  computed  chlorophyll  to  basic  particle 
behaviors.  The  second  examined  sensitivity  of  particle  location  and 
computed  chlorophyll  to  hydrology. 

The  authors  performed  three  experiments  that  initiated  the  investigation 
into  the  potential  effects  of  buoyancy  regulation  on  particle  residence  time 
and  location  and  on  chlorophyll  concentration: 


1.  Particles  rise  at  the  rate  of  1  m  d1. 
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2.  Particles  sink  at  the  rate  of  l  m  cK 

3.  Particles  sink  at  the  rate  of  1  m  d1  in  the  LE  portion  of  the  estuary  and  rise 
at  the  rate  1  m  d1  in  the  RET  portion  of  the  estuary. 

Buoyant  (rising)  particles  were  rapidly  carried  out  of  the  estuary  by  surface 
currents.  Few  particles  left  the  system  in  the  runs  in  which  particles  sink 
and  in  which  particles  rise  and  sink.  Instead  particles  accumulated  in  the 
mid-depth  and  bottom  of  the  RET  reach.  These  runs  indicated  that  sinking 
tends  to  increase  the  residence  time  of  particles  in  the  system. 

Chlorophyll  concentrations  were  examined  by  seasonally  averaging 
computations  along  the  longitudinal  axis  of  the  river.  Particle  buoyancy 
diminished  the  chlorophyll  concentration  along  the  bottom  in  both  winter 
and  spring  relative  to  the  base  run  with  neutrally  buoyant  particles.  The  run 
with  buoyant  particles  provided  superior  surface  chlorophyll  computations 
in  winter  and  computations  equivalent  to  the  neutrallybuoyant  run  in 
spring.  The  run  in  which  particles  rise  and  sink  provided  a  surface 
chlorophyll  peak  90  upstream  from  the  river  mouth  and  partially  validated 
the  conceptual  model  of  algal  blooms  that  prompted  this  study.  Elsewhere, 
however,  results  from  the  run  in  which  algae  rise  and  sink  were 
unsatisfactory. 

The  year  1994,  which  was  selected  for  the  base  runs,  is  a  year  of  average 
hydrology  in  the  Chesapeake  Bay  system  and  of  above-average  runoff  in 
the  Potomac  River.  The  year  1999  is,  by  contrast,  a  dry  year  in  the 
Chesapeake  Bay  system.  Average  flow  at  the  Potomac  fall  line  for  the 
period  November  1998  through  May  1999  is  208  m3  s1,  compared  to  802 
m3  s1  for  the  period  November  1993  through  May  1994.  To  test  the 
sensitivity  of  model  results  to  hydrology,  the  model  of  algae  as  particles 
was  run  for  the  period  November  1998  through  May  1999,  using  the  same 
initial  conditions  and  parameter  set  as  the  base  run. 

Results  from  the  dry  year  showed  more  particles  upstream  of  km  60  than 
the  wet  year,  while  the  wet  year  had  more  particles  in  the  downstream 
portion  of  the  estuary  than  the  dry  year.  The  change  in  spatial  distribution 
of  particles  was  entirely  consistent  with  the  relative  magnitudes  of  the  net 
flows  in  the  two  years.  Chlorophyll  results  from  the  two  years  were  similar. 
Both  runs  demonstrated  similar  periodic  behavior,  similar  spatial  distribu¬ 
tions,  and  similar  comparisons  to  observations.  The  large  alteration  in 
hydrology,  from  1994  to  1999,  induced  no  aberrant  behavior  in  the  model 
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and  confirmed  that  the  model  is  suitable  for  further  investigation  and 
application. 

Conclusions 

The  characteristics  and  applicability  of  the  particle  approach  to  modeling 
algal  blooms  are  now  established.  The  particle  approach  yields  results  that 
differ  distinctly  from  a  conventional  model  of  algae  as  a  dissolved  sub¬ 
stance.  Several  aspects  of  the  particle  approach  are  superior  to  the 
conventional  approach.  In  particular,  the  particle  approach  reproduces 
elevated  chlorophyll  concentrations  observed  at  great  depths  in  the 
Potomac  River.  The  particle  approach  is  not  clearly  superior  to  the  conven¬ 
tional  approach  in  reproducing  the  population  of  chlorophyll  observations. 
Results  from  the  conventional  model  reflect  intense  effort  to  tune  the 
model  to  observations.  For  comparison  purposes,  the  particle  model  used 
an  identical  parameter  set.  Improved  results  may  be  expected  from  the 
particle  model  following  development  of  a  specific  parameter  set.  Now  that 
the  particle  model  has  been  formulated,  explored,  and  documented, 
dedicated  applications  to  blooms  in  the  Potomac  River  and  elsewhere  are 
recommended. 
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